]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implement and test the alternating form for faces.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Thu, 16 Jul 2009 17:54:40 +0000 (17:54 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Thu, 16 Jul 2009 17:54:40 +0000 (17:54 +0000)
git-svn-id: https://svn.dealii.org/trunk@19098 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/source/geometry_info.cc
tests/base/geometry_info_7.cc [new file with mode: 0644]
tests/base/geometry_info_7/cmp/generic [new file with mode: 0644]

index 81e07c7b3c376f00db5667c0b96c48f6dc69bd80..5f61358bc119179661e7cc4e35bae342841b86fd 100644 (file)
@@ -1684,6 +1684,17 @@ namespace internal
     }
 
 
+    inline
+    Tensor<1,3>
+    wedge_product_of_columns (const Tensor<1,3> (&derivative)[2])
+    {
+      Tensor<1,3> result;
+      cross_product (result, derivative[0], derivative[1]);
+
+      return result;
+    }
+    
+
                                     /**
                                      * Alternating form for quads in 3d.
                                      */
@@ -1695,11 +1706,40 @@ namespace internal
     {
                                       // for each of the vertices,
                                       // form the scaled normal
-                                      // vector. since the mapping is
-                                      // linear, all normals are the
-                                      // same
-      (void)vertices;
-      (void)forms;
+                                      // vector. to do so, we need to
+                                      // see how the infinitesimal
+                                      // vectors (d\xi_1,0) and (0,d\xi_2)
+                                      // are transformed into
+                                      // spacedim-dimensional space
+                                      // and then form their cross
+                                      // product. to this end, note that
+                                      //    \vec x = sum_i \vec v_i phi_i(\vec xi)
+                                      // so the transformed vectors are
+                                      //   [x(\xi+(d\xi_1,0))-x(\xi)]/d\xi_1
+                                      // and
+                                      //   [x(\xi+(0,d\xi_2))-x(\xi)]/d\xi_2
+                                      // which boils down to the columns
+                                      // of the 3x2 matrix \grad_\xi x(\xi)
+      const unsigned int dim      = 2;
+      const unsigned int spacedim = 3;
+
+      for (unsigned int i=0; i<dealii::GeometryInfo<dim>::vertices_per_cell; ++i)
+       {
+         Tensor<1,spacedim> derivatives[dim];
+      
+         for (unsigned int j=0; j<dealii::GeometryInfo<dim>::vertices_per_cell; ++j)
+           {
+             const Tensor<1,dim> grad_phi_j
+               = (dealii::GeometryInfo<dim>::
+                  d_linear_shape_function_gradient (dealii::GeometryInfo<dim>::
+                                                    unit_cell_vertex(i),
+                                                    j));
+             for (unsigned int l=0; l<dim; ++l)
+               derivatives[l] += vertices[j] * grad_phi_j[l];
+           }
+
+         forms[i] = wedge_product_of_columns (derivatives);
+       }
     }
   }
 }
diff --git a/tests/base/geometry_info_7.cc b/tests/base/geometry_info_7.cc
new file mode 100644 (file)
index 0000000..d0f8ae5
--- /dev/null
@@ -0,0 +1,200 @@
+//----------------------------  geometry_info_7.cc  ---------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2003, 2004, 2005, 2006, 2009 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.
+//
+//----------------------------  geometry_info_7.cc  ---------------------------
+
+
+// check GeometryInfo::alternating_form_at_vertices for the faces of cells
+
+#include "../tests.h"
+#include <base/logstream.h>
+#include <base/geometry_info.h>
+
+#include <fstream>
+#include <cstdlib>
+
+
+template <int dim>
+void test ()
+{
+  deallog << "Checking in " << dim << "d" << std::endl;
+
+                                  // check the determinant of the
+                                  // transformation for the reference
+                                  // cell. the determinant should be one in
+                                  // that case
+  {
+    Point<dim> vertices[GeometryInfo<dim>::vertices_per_cell];
+    for (unsigned int v=0;v<GeometryInfo<dim>::vertices_per_cell;++v)
+      vertices[v] = GeometryInfo<dim>::unit_cell_vertex(v);
+
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+      {
+       Point<dim> face_vertices[GeometryInfo<dim>::vertices_per_face];
+       for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_face; ++v)
+         face_vertices[v] = vertices[GeometryInfo<dim>::face_to_cell_vertices (f, v)];
+       
+       Tensor<1,dim> alternating_forms[GeometryInfo<dim>::vertices_per_face];
+       GeometryInfo<dim-1>::alternating_form_at_vertices (face_vertices,
+                                                          alternating_forms);
+       for (unsigned int v=0;v<GeometryInfo<dim>::vertices_per_face;++v)
+         {
+           deallog << "Reference cell: face " << f << ": " << alternating_forms[v]
+                   << std::endl;
+           Assert (alternating_forms[v].norm() == 1, ExcInternalError());
+         }
+      }
+  }
+  
+                                  // try the same, but move squash the cell
+                                  // in the x-direction by a factor of 10
+  {
+    Point<dim> vertices[GeometryInfo<dim>::vertices_per_cell];
+    for (unsigned int v=0;v<GeometryInfo<dim>::vertices_per_cell;++v)
+      {
+       vertices[v] = GeometryInfo<dim>::unit_cell_vertex(v);
+       vertices[v][0] /= 10;
+      }
+
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+      {
+       Point<dim> face_vertices[GeometryInfo<dim>::vertices_per_face];
+       for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_face; ++v)
+         face_vertices[v] = vertices[GeometryInfo<dim>::face_to_cell_vertices (f, v)];
+       Tensor<1,dim> alternating_forms[GeometryInfo<dim>::vertices_per_face];
+       GeometryInfo<dim-1>::alternating_form_at_vertices (face_vertices,
+                                                          alternating_forms);
+       for (unsigned int v=0;v<GeometryInfo<dim>::vertices_per_face;++v)
+         {
+           deallog << "Squashed cell: face " << f << ": " << alternating_forms[v]
+                   << std::endl;
+                                            // faces 0,1 should be
+                                            // unaffected, but all
+                                            // other faces are
+                                            // squashed
+           if (f < 2)
+             Assert (alternating_forms[v].norm() == 1, ExcInternalError())
+           else
+             Assert (alternating_forms[v].norm() == 0.1, ExcInternalError());
+         }
+      }
+  }
+  
+
+                                  // try the same, but move squash the cell
+                                  // in the x-direction by a factor of 10 and
+                                  // rotate it around the z-axis (unless in
+                                  // 1d)
+  {
+    Point<dim> vertices[GeometryInfo<dim>::vertices_per_cell];
+    for (unsigned int v=0;v<GeometryInfo<dim>::vertices_per_cell;++v)
+      {
+       vertices[v] = GeometryInfo<dim>::unit_cell_vertex(v);
+       vertices[v][0] /= 10;
+
+       if (dim > 1)
+         {
+           std::swap (vertices[v][0], vertices[v][1]);
+           vertices[v][1] *= -1;
+         }
+      }
+
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+      {
+       Point<dim> face_vertices[GeometryInfo<dim>::vertices_per_face];
+       for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_face; ++v)
+         face_vertices[v] = vertices[GeometryInfo<dim>::face_to_cell_vertices (f, v)];
+       Tensor<1,dim> alternating_forms[GeometryInfo<dim>::vertices_per_face];
+       GeometryInfo<dim-1>::alternating_form_at_vertices (face_vertices,
+                                                          alternating_forms);
+       for (unsigned int v=0;v<GeometryInfo<dim>::vertices_per_face;++v)
+         {
+           deallog << "Squashed+rotated cell: face " << f << ": " << alternating_forms[v]
+                   << std::endl;
+
+                                            // in 2d and 3d, faces
+                                            // 0,1 should be
+                                            // unaffected (just like
+                                            // for the squashed cell,
+                                            // the rotation has
+                                            // nothing to do with
+                                            // face numbers though
+                                            // the direction of the
+                                            // alternating form
+                                            // vector would have
+                                            // rotated along)
+           if (f<2)
+             Assert (alternating_forms[v].norm() == 1, ExcInternalError())
+           else
+             Assert (alternating_forms[v].norm() == 0.1, ExcInternalError());
+         }
+      }
+  }
+  
+                                  // pinched cell
+  {
+    Point<dim> vertices[GeometryInfo<dim>::vertices_per_cell];
+    for (unsigned int v=0;v<GeometryInfo<dim>::vertices_per_cell;++v)
+      vertices[v] = GeometryInfo<dim>::unit_cell_vertex(v);
+    vertices[1] /= 10;
+    
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+      {
+       Point<dim> face_vertices[GeometryInfo<dim>::vertices_per_face];
+       for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_face; ++v)
+         face_vertices[v] = vertices[GeometryInfo<dim>::face_to_cell_vertices (f, v)];
+       Tensor<1,dim> alternating_forms[GeometryInfo<dim>::vertices_per_face];
+       GeometryInfo<dim-1>::alternating_form_at_vertices (face_vertices,
+                                                          alternating_forms);
+       for (unsigned int v=0;v<GeometryInfo<dim>::vertices_per_face;++v)
+         deallog << "Pinched cell: face " << f << ": " << alternating_forms[v]
+                 << std::endl;
+      }
+  }
+  
+
+                                  // inverted cell
+  {
+    Point<dim> vertices[GeometryInfo<dim>::vertices_per_cell];
+    for (unsigned int v=0;v<GeometryInfo<dim>::vertices_per_cell;++v)
+      vertices[v] = GeometryInfo<dim>::unit_cell_vertex(v);
+    std::swap (vertices[0], vertices[1]);
+    
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+      {
+       Point<dim> face_vertices[GeometryInfo<dim>::vertices_per_face];
+       for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_face; ++v)
+         face_vertices[v] = vertices[GeometryInfo<dim>::face_to_cell_vertices (f, v)];
+
+       Tensor<1,dim> alternating_forms[GeometryInfo<dim>::vertices_per_face];
+       GeometryInfo<dim-1>::alternating_form_at_vertices (face_vertices,
+                                                          alternating_forms);
+       for (unsigned int v=0;v<GeometryInfo<dim>::vertices_per_face;++v)
+         deallog << "Inverted cell: face " << f << ": " << alternating_forms[v]
+                 << std::endl;
+      }
+  }
+}
+
+
+
+int main () 
+{
+  std::ofstream logfile("geometry_info_7/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  test<2> ();
+  test<3> ();
+
+  return 0;
+}
diff --git a/tests/base/geometry_info_7/cmp/generic b/tests/base/geometry_info_7/cmp/generic
new file mode 100644 (file)
index 0000000..57f659b
--- /dev/null
@@ -0,0 +1,163 @@
+
+DEAL::Checking in 2d
+DEAL::Reference cell: face 0: 1.00000 0.00000
+DEAL::Reference cell: face 0: 1.00000 0.00000
+DEAL::Reference cell: face 1: 1.00000 0.00000
+DEAL::Reference cell: face 1: 1.00000 0.00000
+DEAL::Reference cell: face 2: 0.00000 -1.00000
+DEAL::Reference cell: face 2: 0.00000 -1.00000
+DEAL::Reference cell: face 3: 0.00000 -1.00000
+DEAL::Reference cell: face 3: 0.00000 -1.00000
+DEAL::Squashed cell: face 0: 1.00000 0.00000
+DEAL::Squashed cell: face 0: 1.00000 0.00000
+DEAL::Squashed cell: face 1: 1.00000 0.00000
+DEAL::Squashed cell: face 1: 1.00000 0.00000
+DEAL::Squashed cell: face 2: 0.00000 -0.100000
+DEAL::Squashed cell: face 2: 0.00000 -0.100000
+DEAL::Squashed cell: face 3: 0.00000 -0.100000
+DEAL::Squashed cell: face 3: 0.00000 -0.100000
+DEAL::Squashed+rotated cell: face 0: 0.00000 -1.00000
+DEAL::Squashed+rotated cell: face 0: 0.00000 -1.00000
+DEAL::Squashed+rotated cell: face 1: 0.00000 -1.00000
+DEAL::Squashed+rotated cell: face 1: 0.00000 -1.00000
+DEAL::Squashed+rotated cell: face 2: -0.100000 0.00000
+DEAL::Squashed+rotated cell: face 2: -0.100000 0.00000
+DEAL::Squashed+rotated cell: face 3: -0.100000 0.00000
+DEAL::Squashed+rotated cell: face 3: -0.100000 0.00000
+DEAL::Pinched cell: face 0: 1.00000 0.00000
+DEAL::Pinched cell: face 0: 1.00000 0.00000
+DEAL::Pinched cell: face 1: 1.00000 -0.900000
+DEAL::Pinched cell: face 1: 1.00000 -0.900000
+DEAL::Pinched cell: face 2: 0.00000 -0.100000
+DEAL::Pinched cell: face 2: 0.00000 -0.100000
+DEAL::Pinched cell: face 3: 0.00000 -1.00000
+DEAL::Pinched cell: face 3: 0.00000 -1.00000
+DEAL::Inverted cell: face 0: 1.00000 1.00000
+DEAL::Inverted cell: face 0: 1.00000 1.00000
+DEAL::Inverted cell: face 1: 1.00000 -1.00000
+DEAL::Inverted cell: face 1: 1.00000 -1.00000
+DEAL::Inverted cell: face 2: 0.00000 1.00000
+DEAL::Inverted cell: face 2: 0.00000 1.00000
+DEAL::Inverted cell: face 3: 0.00000 -1.00000
+DEAL::Inverted cell: face 3: 0.00000 -1.00000
+DEAL::Checking in 3d
+DEAL::Reference cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Reference cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Reference cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Reference cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Reference cell: face 1: 1.00000 0.00000 0.00000
+DEAL::Reference cell: face 1: 1.00000 0.00000 0.00000
+DEAL::Reference cell: face 1: 1.00000 0.00000 0.00000
+DEAL::Reference cell: face 1: 1.00000 0.00000 0.00000
+DEAL::Reference cell: face 2: 0.00000 1.00000 0.00000
+DEAL::Reference cell: face 2: 0.00000 1.00000 0.00000
+DEAL::Reference cell: face 2: 0.00000 1.00000 0.00000
+DEAL::Reference cell: face 2: 0.00000 1.00000 0.00000
+DEAL::Reference cell: face 3: 0.00000 1.00000 0.00000
+DEAL::Reference cell: face 3: 0.00000 1.00000 0.00000
+DEAL::Reference cell: face 3: 0.00000 1.00000 0.00000
+DEAL::Reference cell: face 3: 0.00000 1.00000 0.00000
+DEAL::Reference cell: face 4: 0.00000 0.00000 1.00000
+DEAL::Reference cell: face 4: 0.00000 0.00000 1.00000
+DEAL::Reference cell: face 4: 0.00000 0.00000 1.00000
+DEAL::Reference cell: face 4: 0.00000 0.00000 1.00000
+DEAL::Reference cell: face 5: 0.00000 0.00000 1.00000
+DEAL::Reference cell: face 5: 0.00000 0.00000 1.00000
+DEAL::Reference cell: face 5: 0.00000 0.00000 1.00000
+DEAL::Reference cell: face 5: 0.00000 0.00000 1.00000
+DEAL::Squashed cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Squashed cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Squashed cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Squashed cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Squashed cell: face 1: 1.00000 0.00000 0.00000
+DEAL::Squashed cell: face 1: 1.00000 0.00000 0.00000
+DEAL::Squashed cell: face 1: 1.00000 0.00000 0.00000
+DEAL::Squashed cell: face 1: 1.00000 0.00000 0.00000
+DEAL::Squashed cell: face 2: 0.00000 0.100000 0.00000
+DEAL::Squashed cell: face 2: 0.00000 0.100000 0.00000
+DEAL::Squashed cell: face 2: 0.00000 0.100000 0.00000
+DEAL::Squashed cell: face 2: 0.00000 0.100000 0.00000
+DEAL::Squashed cell: face 3: 0.00000 0.100000 0.00000
+DEAL::Squashed cell: face 3: 0.00000 0.100000 0.00000
+DEAL::Squashed cell: face 3: 0.00000 0.100000 0.00000
+DEAL::Squashed cell: face 3: 0.00000 0.100000 0.00000
+DEAL::Squashed cell: face 4: 0.00000 0.00000 0.100000
+DEAL::Squashed cell: face 4: 0.00000 0.00000 0.100000
+DEAL::Squashed cell: face 4: 0.00000 0.00000 0.100000
+DEAL::Squashed cell: face 4: 0.00000 0.00000 0.100000
+DEAL::Squashed cell: face 5: 0.00000 0.00000 0.100000
+DEAL::Squashed cell: face 5: 0.00000 0.00000 0.100000
+DEAL::Squashed cell: face 5: 0.00000 0.00000 0.100000
+DEAL::Squashed cell: face 5: 0.00000 0.00000 0.100000
+DEAL::Squashed+rotated cell: face 0: 0.00000 -1.00000 0.00000
+DEAL::Squashed+rotated cell: face 0: 0.00000 -1.00000 0.00000
+DEAL::Squashed+rotated cell: face 0: 0.00000 -1.00000 0.00000
+DEAL::Squashed+rotated cell: face 0: 0.00000 -1.00000 0.00000
+DEAL::Squashed+rotated cell: face 1: 0.00000 -1.00000 0.00000
+DEAL::Squashed+rotated cell: face 1: 0.00000 -1.00000 0.00000
+DEAL::Squashed+rotated cell: face 1: 0.00000 -1.00000 0.00000
+DEAL::Squashed+rotated cell: face 1: 0.00000 -1.00000 0.00000
+DEAL::Squashed+rotated cell: face 2: 0.100000 0.00000 0.00000
+DEAL::Squashed+rotated cell: face 2: 0.100000 0.00000 0.00000
+DEAL::Squashed+rotated cell: face 2: 0.100000 0.00000 0.00000
+DEAL::Squashed+rotated cell: face 2: 0.100000 0.00000 0.00000
+DEAL::Squashed+rotated cell: face 3: 0.100000 0.00000 0.00000
+DEAL::Squashed+rotated cell: face 3: 0.100000 0.00000 0.00000
+DEAL::Squashed+rotated cell: face 3: 0.100000 0.00000 0.00000
+DEAL::Squashed+rotated cell: face 3: 0.100000 0.00000 0.00000
+DEAL::Squashed+rotated cell: face 4: 0.00000 0.00000 0.100000
+DEAL::Squashed+rotated cell: face 4: 0.00000 0.00000 0.100000
+DEAL::Squashed+rotated cell: face 4: 0.00000 0.00000 0.100000
+DEAL::Squashed+rotated cell: face 4: 0.00000 0.00000 0.100000
+DEAL::Squashed+rotated cell: face 5: 0.00000 0.00000 0.100000
+DEAL::Squashed+rotated cell: face 5: 0.00000 0.00000 0.100000
+DEAL::Squashed+rotated cell: face 5: 0.00000 0.00000 0.100000
+DEAL::Squashed+rotated cell: face 5: 0.00000 0.00000 0.100000
+DEAL::Pinched cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Pinched cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Pinched cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Pinched cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Pinched cell: face 1: 1.00000 -0.900000 -0.900000
+DEAL::Pinched cell: face 1: 1.00000 -0.900000 0.00000
+DEAL::Pinched cell: face 1: 1.00000 0.00000 -0.900000
+DEAL::Pinched cell: face 1: 1.00000 0.00000 0.00000
+DEAL::Pinched cell: face 2: 0.00000 0.100000 0.00000
+DEAL::Pinched cell: face 2: 0.00000 1.00000 0.00000
+DEAL::Pinched cell: face 2: 0.00000 0.100000 0.00000
+DEAL::Pinched cell: face 2: 0.00000 1.00000 0.00000
+DEAL::Pinched cell: face 3: 0.00000 1.00000 0.00000
+DEAL::Pinched cell: face 3: 0.00000 1.00000 0.00000
+DEAL::Pinched cell: face 3: 0.00000 1.00000 0.00000
+DEAL::Pinched cell: face 3: 0.00000 1.00000 0.00000
+DEAL::Pinched cell: face 4: 0.00000 0.00000 0.100000
+DEAL::Pinched cell: face 4: 0.00000 0.00000 0.100000
+DEAL::Pinched cell: face 4: 0.00000 0.00000 1.00000
+DEAL::Pinched cell: face 4: 0.00000 0.00000 1.00000
+DEAL::Pinched cell: face 5: 0.00000 0.00000 1.00000
+DEAL::Pinched cell: face 5: 0.00000 0.00000 1.00000
+DEAL::Pinched cell: face 5: 0.00000 0.00000 1.00000
+DEAL::Pinched cell: face 5: 0.00000 0.00000 1.00000
+DEAL::Inverted cell: face 0: 1.00000 1.00000 1.00000
+DEAL::Inverted cell: face 0: 1.00000 1.00000 0.00000
+DEAL::Inverted cell: face 0: 1.00000 0.00000 1.00000
+DEAL::Inverted cell: face 0: 1.00000 0.00000 0.00000
+DEAL::Inverted cell: face 1: 1.00000 -1.00000 -1.00000
+DEAL::Inverted cell: face 1: 1.00000 -1.00000 0.00000
+DEAL::Inverted cell: face 1: 1.00000 0.00000 -1.00000
+DEAL::Inverted cell: face 1: 1.00000 0.00000 0.00000
+DEAL::Inverted cell: face 2: 0.00000 -1.00000 0.00000
+DEAL::Inverted cell: face 2: 0.00000 1.00000 0.00000
+DEAL::Inverted cell: face 2: 0.00000 -1.00000 0.00000
+DEAL::Inverted cell: face 2: 0.00000 1.00000 0.00000
+DEAL::Inverted cell: face 3: 0.00000 1.00000 0.00000
+DEAL::Inverted cell: face 3: 0.00000 1.00000 0.00000
+DEAL::Inverted cell: face 3: 0.00000 1.00000 0.00000
+DEAL::Inverted cell: face 3: 0.00000 1.00000 0.00000
+DEAL::Inverted cell: face 4: 0.00000 0.00000 -1.00000
+DEAL::Inverted cell: face 4: 0.00000 0.00000 -1.00000
+DEAL::Inverted cell: face 4: 0.00000 0.00000 1.00000
+DEAL::Inverted cell: face 4: 0.00000 0.00000 1.00000
+DEAL::Inverted cell: face 5: 0.00000 0.00000 1.00000
+DEAL::Inverted cell: face 5: 0.00000 0.00000 1.00000
+DEAL::Inverted cell: face 5: 0.00000 0.00000 1.00000
+DEAL::Inverted cell: face 5: 0.00000 0.00000 1.00000

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.