]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Work around a problem of Suns Forte compiler which cant handle the case that an enum...
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 20 Mar 2002 12:46:21 +0000 (12:46 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 20 Mar 2002 12:46:21 +0000 (12:46 +0000)
/* ----------------------------------------------- */
/* Problem 10 -- Members with the same name as the */
/* enclosing namespace.                            */
namespace NS3 {
  class NS3 {};
};
NS3::NS3 ns3;

Thus rename the IteratorState enum in the IteratorState namespace into IteratorState_s_.

git-svn-id: https://svn.dealii.org/trunk@5599 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/tria_accessor.h
deal.II/deal.II/include/grid/tria_accessor.templates.h
deal.II/deal.II/include/grid/tria_iterator.h
deal.II/deal.II/include/grid/tria_iterator_base.h
deal.II/examples/step-1/Makefile
deal.II/examples/step-7/Makefile
deal.II/examples/step-7/step-7.cc

index cb6b8789f5208ecde929825ceef6e84521252657..232a1895a26700ccfff140f453bf14c9b233b02b 100644 (file)
@@ -168,7 +168,7 @@ class TriaAccessor
                                      *  can be in, refer to the
                                      *  @ref{TriaRawIterator} documentation.
                                      */
-    IteratorState::IteratorState state () const;
+    IteratorState::IteratorStates state () const;
 
                                     /**
                                      * Return a pointer to the triangulation
index 00a32ba2d58116eed3aeac4753b024e19c256176..e5d17ed87d85dd490ff21ff919395a20c646a894 100644 (file)
@@ -86,7 +86,7 @@ TriaAccessor<dim>::index () const
 
 template <int dim>
 inline
-IteratorState::IteratorState
+IteratorState::IteratorStates
 TriaAccessor<dim>::state () const
 {
   if ((present_level>=0) && (present_index>=0))
index 1bceacd45e173f4c06652b1e712e2c74316daae3..5bb60880f01fa6f473672fb9d242f8fdc0a6312a 100644 (file)
@@ -288,7 +288,7 @@ class TriaRawIterator :
     TriaRawIterator (Triangulation<dim> *parent,
                     const int           level,
                     const int           index,
-                    const typename Accessor::AccessorData *local_data = 0);
+                    const typename AccessorType::AccessorData *local_data = 0);
 
                                     /**
                                      * This is a conversion operator
@@ -465,7 +465,7 @@ class TriaRawIterator :
                                     /**
                                      *  Return the state of the iterator.
                                      */
-    IteratorState::IteratorState state () const;
+    IteratorState::IteratorStates state () const;
 
                                     /**
                                      * Print the iterator to @p{out}. The
@@ -901,7 +901,7 @@ TriaRawIterator<dim,Accessor>::operator -> ()
 
 template <int dim, typename Accessor>
 inline
-IteratorState::IteratorState
+IteratorState::IteratorStates
 TriaRawIterator<dim,Accessor>::state () const
 {
   return accessor.state ();
index 11f3eb5dba0ade29db5e54cabe716561169617a7..7da1415a3dee57c1d3bac3772db69593a087d8b8 100644 (file)
@@ -28,7 +28,7 @@ namespace IteratorState
  *   The three states an iterator can be in: valid, past-the-end and
  *   invalid.
  */
-  enum IteratorState
+  enum IteratorStates
   {
        valid,
        past_the_end,
index 9761232d00e38fa52a00183cec56a6573f39404c..5d76059226ddd16963a73c90b337f777433256bd 100644 (file)
@@ -14,7 +14,7 @@ target = $(basename $(shell echo step-*.cc))
 # run-time checking of parameters and internal states is performed, so
 # you should set this value to `on' while you develop your program,
 # and to `off' when running production computations.
-debug-mode = on
+debug-mode = off
 
 
 # As third field, we need to give the path to the top-level deal.II
index 6360bcca20162bd174d67d4692443f3fc47b97f3..da9acbbdb76972dc6f735b98d969865f886ef03e 100644 (file)
@@ -14,7 +14,7 @@ target = $(basename $(shell echo step-*.cc))
 # run-time checking of parameters and internal states is performed, so
 # you should set this value to `on' while you develop your program,
 # and to `off' when running production computations.
-debug-mode = on
+debug-mode = off
 
 
 # As third field, we need to give the path to the top-level deal.II
index a27f10c6bd59ead3da994e29985817a0db8cb5cb..bc0612742681edea3a4791b6ccfe06f28fa9186f 100644 (file)
@@ -35,6 +35,7 @@
 #include <dofs/dof_accessor.h>
 #include <dofs/dof_tools.h>
 #include <fe/fe_q.h>
+#include <fe/mapping_q.h>
 #include <numerics/matrices.h>
 #include <numerics/error_estimator.h>
 #include <numerics/data_out.h>
@@ -166,6 +167,14 @@ SolutionBase<2>::source_centers[SolutionBase<2>::n_source_centers]
     Point<2>(-0.5, -0.5), 
     Point<2>(+0.5, -0.5)   };
 
+                                // ...and finally for 3d:
+template <>
+const Point<3>
+SolutionBase<3>::source_centers[SolutionBase<3>::n_source_centers]
+= { Point<3>(-0.5, +0.5, -0.5), 
+    Point<3>(-0.5, -0.5, +0.5), 
+    Point<3>(+0.5, -0.5, 0.0)   };
+
                                 // There remains to assign a value to
                                 // the half-width of the
                                 // exponentials. We would like to use
@@ -780,7 +789,8 @@ void LaplaceProblem<dim>::assemble_system ()
                                   // cell (rather than on the unit
                                   // cell) to evaluate the right hand
                                   // side function.
-  FEValues<dim>  fe_values (*fe, quadrature_formula, 
+  MappingQ<dim> mapping(3);
+  FEValues<dim>  fe_values (mapping, *fe, quadrature_formula, 
                            UpdateFlags(update_values    |
                                        update_gradients |
                                        update_q_points  |
@@ -1051,15 +1061,21 @@ void LaplaceProblem<dim>::assemble_system ()
 template <int dim>
 void LaplaceProblem<dim>::solve () 
 {
-  SolverControl           solver_control (1000, 1e-12);
+  SolverControl           solver_control (30, 1e-4);
   PrimitiveVectorMemory<> vector_memory;
   SolverCG<>              cg (solver_control, vector_memory);
 
   PreconditionSSOR<> preconditioner;
   preconditioner.initialize(system_matrix, 1.2);
 
-  cg.solve (system_matrix, solution, system_rhs,
-           preconditioner);
+  try
+    {
+      cg.solve (system_matrix, solution, system_rhs,
+               preconditioner);
+    }
+  catch (...)
+    {
+    };
 
   hanging_node_constraints.distribute (solution);
 };
@@ -1322,7 +1338,11 @@ void LaplaceProblem<dim>::process_solution (const unsigned int cycle)
 template <int dim>
 void LaplaceProblem<dim>::run () 
 {
-  for (unsigned int cycle=0; cycle<7; ++cycle)
+  const unsigned int n_cycles = (refinement_mode == global_refinement ?
+                                (dim == 2 ? 7 : 3)
+                                :
+                                (dim == 2 ? 12 : 7));
+  for (unsigned int cycle=0; cycle<n_cycles; ++cycle)
     {
                                       // The first action in each
                                       // iteration of the outer loop

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.