<commit_before>//-----------------------------------------------------------------
//    $Id$
//    Version: $Name$
//
//    Copyright (C) 1998, 1999, 2000, 2001 by the deal.II authors
//
//    This file is subject to QPL and may not be  distributed
//    without copyright and license information. Please refer
//    to the file deal.II/doc/license.html for the  text  and
//    further information on this license.
//
//-----------------------------------------------------------------

#include <base/tensor.h>
#include <base/quadrature.h>
#include <lac/full_matrix.h>
#include <grid/tria.h>
#include <grid/tria_iterator.h>
#include <dofs/dof_accessor.h>
#include <fe/mapping_cartesian.h>
#include <fe/fe_values.h>

#include <cmath>
#include <algorithm>




template <int dim>
const unsigned int MappingCartesian<dim>::invalid_face_number;



template<int dim>
MappingCartesian<dim>::InternalData::InternalData (const Quadrature<dim>& q)
		:
		quadrature_points (q.get_points ())
{}



template <int dim>
UpdateFlags
MappingCartesian<dim>::update_once (const UpdateFlags) const
{
  return update_default;
}



template <int dim>
UpdateFlags
MappingCartesian<dim>::update_each (const UpdateFlags in) const
{
  UpdateFlags out = in;
  if (out & update_boundary_forms)
    out |= update_normal_vectors;
  
  return out;
}



template <int dim>
typename Mapping<dim>::InternalDataBase *
MappingCartesian<dim>::get_data (const UpdateFlags      flags,
				 const Quadrature<dim> &q) const
{
  Assert (flags & update_normal_vectors == 0, ExcNotImplemented());
  InternalData* data = new InternalData (q);
  return data;
}



template <int dim>
typename Mapping<dim>::InternalDataBase *
MappingCartesian<dim>::get_face_data (const UpdateFlags,
				      const Quadrature<dim-1>& quadrature) const
{
  QProjector<dim> q (quadrature, false);
  InternalData* data = new InternalData (q);
  return data;
}



template <int dim>
typename Mapping<dim>::InternalDataBase *
MappingCartesian<dim>::get_subface_data (const UpdateFlags,
					 const Quadrature<dim-1> &quadrature) const
{
  QProjector<dim> q (quadrature, true);
  InternalData* data = new InternalData (q);
  return data;
}




template <int dim>
void
MappingCartesian<dim>::compute_fill (const typename DoFHandler<dim>::cell_iterator &cell,
				     const unsigned int        face_no,
				     const unsigned int        sub_no,
				     InternalData             &data,
				     std::vector<Point<dim> > &quadrature_points,
				     std::vector<Point<dim> > &normal_vectors) const
{
//TODO:[GK] does it make sense to query the update flags in data if we did not set any in the get_*_data functions?
  UpdateFlags update_flags(data.current_update_flags());

  const unsigned int npts = quadrature_points.size ();
  unsigned int offset = 0;
  
  if (face_no != invalid_face_number)
    {
				       // Add 1 on both sides of
				       // assertion to avoid compiler
				       // warning about testing
				       // unsigned int < 0 in 1d.
      Assert (face_no+1 < GeometryInfo<dim>::faces_per_cell+1,
	      ExcIndexRange (face_no, 0, GeometryInfo<dim>::faces_per_cell));

      if (sub_no == invalid_face_number)
					 // called from FEFaceValues
	offset = face_no * quadrature_points.size();
      else
	{
					   // called from FESubfaceValue
	  Assert (sub_no+1 < GeometryInfo<dim>::subfaces_per_face+1,
		  ExcIndexRange (sub_no, 0, GeometryInfo<dim>::subfaces_per_face));
	  offset = (face_no * GeometryInfo<dim>::subfaces_per_face + sub_no)
		   * quadrature_points.size();
	}
    }
  else
				     // invalid face number, so
				     // subface should be invalid as
				     // well
    Assert (sub_no == invalid_face_number, ExcInternalError());

				   // let @p{start} be the origin of a
				   // local coordinate system. it is
				   // chosen as the (lower) left
				   // vertex
  const Point<dim> start = cell->vertex(0);
  
				   // Compute start point and sizes
				   // along axes.  Strange vertex
				   // numbering makes this complicated
				   // again.
  switch (dim)
    {
      case 1:
	data.length[0] = cell->vertex(1)(0) - start(0);
	break;
      case 2:
	data.length[0] = cell->vertex(1)(0) - start(0);
	data.length[1] = cell->vertex(3)(1) - start(1);
	break;
      case 3:
	data.length[0] = cell->vertex(1)(0) - start(0);
	data.length[1] = cell->vertex(4)(1) - start(1);
	data.length[2] = cell->vertex(3)(2) - start(2);
	break;
      default:
	Assert(false, ExcNotImplemented());
    }
  

				   // transform quadrature point. this
				   // is obtained simply by scaling
				   // unit coordinates with lengths in
				   // each direction
  if (update_flags & update_q_points)
    {
      Assert (quadrature_points.size() == npts,
	      ExcDimensionMismatch(quadrature_points.size(), npts));
      for (unsigned int i=0; i<npts; ++i)
	{
	  quadrature_points[i] = start;
	  for (unsigned int d=0; d<dim; ++d)
	    quadrature_points[i](d) += data.length[d] *
				       data.quadrature_points[i+offset](d);
	}
    }


				   // compute normal vectors. since
				   // cells are aligned to coordinate
				   // axes, they are simply vectors
				   // with exactly one entry equal to
				   // 1 or -1
  if (update_flags & update_normal_vectors)
    {
      Assert (normal_vectors.size() == npts,
	      ExcDimensionMismatch(normal_vectors.size(), npts));
      Point<dim> n;
      switch (100*dim+face_no)
	{
					   // 2D
	  case 200:
		n (1) = -1.;
		break;
	  case 201:
		n (0) = 1.;
		break;
	  case 202:
		n (1) = 1.;
		break;
	  case 203:
		n (0) = 1.;
		break;
					   // 3D<commit_msg>bug in normal vectors fixed<commit_after>//-----------------------------------------------------------------
//    $Id$
//    Version: $Name$
//
//    Copyright (C) 1998, 1999, 2000, 2001 by the deal.II authors
//
//    This file is subject to QPL and may not be  distributed
//    without copyright and license information. Please refer
//    to the file deal.II/doc/license.html for the  text  and
//    further information on this license.
//
//-----------------------------------------------------------------

#include <base/tensor.h>
#include <base/quadrature.h>
#include <lac/full_matrix.h>
#include <grid/tria.h>
#include <grid/tria_iterator.h>
#include <dofs/dof_accessor.h>
#include <fe/mapping_cartesian.h>
#include <fe/fe_values.h>

#include <cmath>
#include <algorithm>




template <int dim>
const unsigned int MappingCartesian<dim>::invalid_face_number;



template<int dim>
MappingCartesian<dim>::InternalData::InternalData (const Quadrature<dim>& q)
		:
		quadrature_points (q.get_points ())
{}



template <int dim>
UpdateFlags
MappingCartesian<dim>::update_once (const UpdateFlags) const
{
  return update_default;
}



template <int dim>
UpdateFlags
MappingCartesian<dim>::update_each (const UpdateFlags in) const
{
  UpdateFlags out = in;
  if (out & update_boundary_forms)
    out |= update_normal_vectors;
  
  return out;
}



template <int dim>
typename Mapping<dim>::InternalDataBase *
MappingCartesian<dim>::get_data (const UpdateFlags      update_flags,
				 const Quadrature<dim> &q) const
{
  //  Assert (flags & update_normal_vectors == 0, ExcNotImplemented());

  InternalData* data = new InternalData (q);

  data->update_once = update_once(update_flags);
  data->update_each = update_each(update_flags);
  data->update_flags = data->update_once | data->update_each;

  return data;
}



template <int dim>
typename Mapping<dim>::InternalDataBase *
MappingCartesian<dim>::get_face_data (const UpdateFlags update_flags,
				      const Quadrature<dim-1>& quadrature) const
{
  QProjector<dim> q (quadrature, false);
  InternalData* data = new InternalData (q);

  data->update_once = update_once(update_flags);
  data->update_each = update_each(update_flags);
  data->update_flags = data->update_once | data->update_each;

  return data;
}



template <int dim>
typename Mapping<dim>::InternalDataBase *
MappingCartesian<dim>::get_subface_data (const UpdateFlags update_flags,
					 const Quadrature<dim-1> &quadrature) const
{
  QProjector<dim> q (quadrature, true);
  InternalData* data = new InternalData (q);

  data->update_once = update_once(update_flags);
  data->update_each = update_each(update_flags);
  data->update_flags = data->update_once | data->update_each;

  return data;
}




template <int dim>
void
MappingCartesian<dim>::compute_fill (const typename DoFHandler<dim>::cell_iterator &cell,
				     const unsigned int        face_no,
				     const unsigned int        sub_no,
				     InternalData             &data,
				     std::vector<Point<dim> > &quadrature_points,
				     std::vector<Point<dim> > &normal_vectors) const
{
  UpdateFlags update_flags(data.current_update_flags());

  const unsigned int npts = quadrature_points.size ();
  unsigned int offset = 0;
  
  if (face_no != invalid_face_number)
    {
				       // Add 1 on both sides of
				       // assertion to avoid compiler
				       // warning about testing
				       // unsigned int < 0 in 1d.
      Assert (face_no+1 < GeometryInfo<dim>::faces_per_cell+1,
	      ExcIndexRange (face_no, 0, GeometryInfo<dim>::faces_per_cell));

      if (sub_no == invalid_face_number)
					 // called from FEFaceValues
	offset = face_no * quadrature_points.size();
      else
	{
					   // called from FESubfaceValue
	  Assert (sub_no+1 < GeometryInfo<dim>::subfaces_per_face+1,
		  ExcIndexRange (sub_no, 0, GeometryInfo<dim>::subfaces_per_face));
	  offset = (face_no * GeometryInfo<dim>::subfaces_per_face + sub_no)
		   * quadrature_points.size();
	}
    }
  else
				     // invalid face number, so
				     // subface should be invalid as
				     // well
    Assert (sub_no == invalid_face_number, ExcInternalError());

				   // let @p{start} be the origin of a
				   // local coordinate system. it is
				   // chosen as the (lower) left
				   // vertex
  const Point<dim> start = cell->vertex(0);
  
				   // Compute start point and sizes
				   // along axes.  Strange vertex
				   // numbering makes this complicated
				   // again.
  switch (dim)
    {
      case 1:
	data.length[0] = cell->vertex(1)(0) - start(0);
	break;
      case 2:
	data.length[0] = cell->vertex(1)(0) - start(0);
	data.length[1] = cell->vertex(3)(1) - start(1);
	break;
      case 3:
	data.length[0] = cell->vertex(1)(0) - start(0);
	data.length[1] = cell->vertex(4)(1) - start(1);
	data.length[2] = cell->vertex(3)(2) - start(2);
	break;
      default:
	Assert(false, ExcNotImplemented());
    }
  

				   // transform quadrature point. this
				   // is obtained simply by scaling
				   // unit coordinates with lengths in
				   // each direction
  if (update_flags & update_q_points)
    {
      Assert (quadrature_points.size() == npts,
	      ExcDimensionMismatch(quadrature_points.size(), npts));
      for (unsigned int i=0; i<npts; ++i)
	{
	  quadrature_points[i] = start;
	  for (unsigned int d=0; d<dim; ++d)
	    quadrature_points[i](d) += data.length[d] *
				       data.quadrature_points[i+offset](d);
	}
    }


				   // compute normal vectors. since
				   // cells are aligned to coordinate
				   // axes, they are simply vectors
				   // with exactly one entry equal to
				   // 1 or -1
  if (update_flags & update_normal_vectors)
    {
      Assert (normal_vectors.size() == npts,
	      ExcDimensionMismatch(normal_vectors.size(), npts));
      Point<dim> n;
      switch (100*dim+face_no)
	{
					   // 2D
	  case 200:
		n (1) = -1.;
		break;
	  case 201:
		n (0) = 1.;
		break;
	  case 202:
		n (1) = 1.;
		break;
	  case 203:
		n (0) = -1.;
		break;
					   // 3D
{'fraction_non_alphanumeric': 0.14384449244060477, 'fraction_numerical': 0.014859611231101511, 'mean_word_length': 2.742644681538959, 'pattern_counts': {'":': 0, '<': 107, '<?xml version=': 0, '>': 126, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '6072627', 'n_tokens_mistral': 4216, 'n_tokens_neox': 3684, 'n_words': 954}