]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Added VectorTools::create_point_source_vector
authorrschulz <rschulz@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 23 May 2006 13:38:59 +0000 (13:38 +0000)
committerrschulz <rschulz@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 23 May 2006 13:38:59 +0000 (13:38 +0000)
git-svn-id: https://svn.dealii.org/trunk@13129 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/Makefile
deal.II/deal.II/include/numerics/vectors.h
deal.II/deal.II/include/numerics/vectors.templates.h
deal.II/deal.II/source/numerics/vectors.cc
deal.II/doc/news/changes.html

index 21a5fb8d11989c08e6c8144ffeeeede391670666..7edf880552f6b8f20d017916193caad230b3f16d 100644 (file)
@@ -125,6 +125,7 @@ $(LIBDIR)/libdeal_II_3d$(static-lib-suffix): $(o-files-3d)
 
 $(LIBDIR)/libdeal_II_1d.g$(shared-lib-suffix): $(go-files-1d)
        @echo "=====deal.II====1d====debug======$(MT)== Linking library:   $(@F)"
+       @echo $(SHLIBLD) $(LDFLAGS) $(SHLIBFLAGS) -o $@ $(go-files-1d) $(deplibs.g)
        @$(SHLIBLD) $(LDFLAGS) $(SHLIBFLAGS) -o $@ $(go-files-1d) $(deplibs.g)
 
 $(LIBDIR)/libdeal_II_1d$(shared-lib-suffix): $(o-files-1d)
@@ -143,10 +144,12 @@ $(LIBDIR)/libdeal_II_2d$(shared-lib-suffix): $(o-files-2d)
 
 $(LIBDIR)/libdeal_II_3d.g$(shared-lib-suffix): $(go-files-3d)
        @echo "=====deal.II====3d====debug======$(MT)== Linking library:   $(@F)"
+       @echo $(SHLIBLD) $(LDFLAGS) $(SHLIBFLAGS) -o $@ $(go-files-3d) $(deplibs.g)
        @$(SHLIBLD) $(LDFLAGS) $(SHLIBFLAGS) -o $@ $(go-files-3d) $(deplibs.g)
 
 $(LIBDIR)/libdeal_II_3d$(shared-lib-suffix): $(o-files-3d)
        @echo "=====deal.II====3d====optimized==$(MT)== Linking library:   $(@F)"
+       @echo $(SHLIBLD) $(LDFLAGS) $(SHLIBFLAGS) -o $@ $(o-files-3d) $(deplibs.o)
        @$(SHLIBLD) $(LDFLAGS) $(SHLIBFLAGS) -o $@ $(o-files-3d) $(deplibs.o)
 
 
index 35dd13ae761c0564497dbf6e68be77b8733666db..3bdfccb08ae4dac2a9c25fbbd07bfbf3cf5cd30f 100644 (file)
@@ -165,6 +165,10 @@ class ConstraintMatrix;
  *   <tt>MatrixCreator::create_*</tt> functions which take a right hand side do,
  *   but without assembling a matrix.
  *
+ * <li> Creation of right hand side vectors for point sources:
+ *   The @p create_point_source_vector function computes the vector
+ *   $f_i = \int_\Omega \delta_0(x-x_0) \phi_i(x) dx$. 
+ *
  * <li> Creation of boundary right hand side vectors: The
  *   @p create_boundary_right_hand_side function computes the vector
  *   $f_i = \int_{\partial\Omega} g(x) \phi_i(x) dx$. This is the
@@ -568,6 +572,32 @@ class VectorTools
                                        Vector<double>        &rhs_vector);
 
                                     /**
+                                     * Create a right hand side
+                                     * vector for a point source at point @p p.
+                                      * Prior content of the
+                                     * given @p rhs_vector vector is
+                                     * deleted.
+                                     *
+                                     * See the general documentation of this
+                                     * class for further information.
+                                     */
+    template <int dim>
+    static void create_point_source_vector(const Mapping<dim>    &mapping,
+                                           const DoFHandler<dim> &dof,
+                                           const Point<dim>      &p,
+                                           Vector<double>        &rhs_vector);
+
+                                    /**
+                                     * Calls the @p create_point_source_vector
+                                     * function, see above, with
+                                     * <tt>mapping=MappingQ1@<dim@>()</tt>.
+                                     */
+    template <int dim>
+    static void create_point_source_vector(const DoFHandler<dim> &dof,
+                                           const Point<dim>      &p,
+                                           Vector<double>        &rhs_vector);
+
+                                     /**
                                      * Create a right hand side
                                      * vector from boundary
                                      * forces. Prior content of the
index 7c665cc40acdbd5c004be0367f066990d2177bec..67068eb5b894157dcc35769399089bfe20b2a422 100644 (file)
@@ -572,12 +572,48 @@ void VectorTools::create_right_hand_side (const DoFHandler<dim>    &dof_handler,
                                          Vector<double>           &rhs_vector)
 {
   Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));
-  static const MappingQ1<dim> mapping;
-  create_right_hand_side(mapping, dof_handler, quadrature,
+  create_right_hand_side(StaticMappingQ1<dim>::mapping, dof_handler, quadrature,
                         rhs_function, rhs_vector);
 }
 
 
+template <int dim>
+void VectorTools::create_point_source_vector (const Mapping<dim>       &mapping,
+                                              const DoFHandler<dim>    &dof_handler,
+                                              const Point<dim>         &p,
+                                              Vector<double>           &rhs_vector)
+{
+   Assert (rhs_vector.size() == dof_handler.n_dofs(),
+           ExcDimensionMismatch(rhs_vector.size(), dof_handler.n_dofs()));
+   rhs_vector = 0;
+
+   std::pair<typename DoFHandler<dim>::active_cell_iterator, Point<dim> >
+      cell_point =
+      GridTools::find_active_cell_around_point (mapping, dof_handler, p);
+
+   Quadrature<dim> q(GeometryInfo<dim>::project_to_unit_cell(cell_point.second));
+
+   FEValues<dim> fe_values(dof_handler.get_fe(), q, UpdateFlags(update_values));
+   fe_values.reinit(cell_point.first);
+
+   const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell;
+
+   std::vector<unsigned int> local_dof_indices (dofs_per_cell);
+   cell_point.first->get_dof_indices (local_dof_indices);
+
+   for(unsigned int i=0; i<dofs_per_cell; i++)
+      rhs_vector(local_dof_indices[i]) =  fe_values.shape_value(i,0);
+}
+
+template <int dim>
+void VectorTools::create_point_source_vector (const DoFHandler<dim>    &dof_handler,
+                                              const Point<dim>         &p,
+                                              Vector<double>           &rhs_vector)
+{
+  Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));
+  create_point_source_vector(StaticMappingQ1<dim>::mapping, dof_handler,
+                             p, rhs_vector);
+}
 
 #if deal_II_dimension != 1
 
index 380eeefb46647f0c0d25c320607b24a40eeb316a..d5d5e6f3a3a22d029d5cf64886600e3286bd1427 100644 (file)
@@ -102,6 +102,17 @@ void VectorTools::create_right_hand_side<deal_II_dimension>
  const Quadrature<deal_II_dimension> &,
  const Function<deal_II_dimension>   &,
  Vector<double>                      &);
+template
+void VectorTools::create_point_source_vector<deal_II_dimension>
+(const Mapping<deal_II_dimension>    &,
+ const DoFHandler<deal_II_dimension> &,
+ const Point<deal_II_dimension>      &,
+ Vector<double>                      &);
+template
+void VectorTools::create_point_source_vector<deal_II_dimension>
+(const DoFHandler<deal_II_dimension> &,
+ const Point<deal_II_dimension>      &,
+ Vector<double>                      &);
 
 #if deal_II_dimension != 1
 template
index d9fc3b0109098a9055999626bfa9a55b2a55799d..eb946c63621fdf2a59e3a33991c2a5b5332118d4 100644 (file)
@@ -646,6 +646,12 @@ inconvenience this causes.
 <h3>deal.II</h3>
 
 <ol>
+  <li> <p>
+       New: Function <code class="class">VectorTools</code>::<code
+       class="class">create_point_source_vector</code> to calculate the projection
+       of a Dirac pulse on the base functions. This models a point source as
+       used in the calculations of Green's functions and can also be used to
+       extract the value of a finite element solution in a single point.
   <li> <p>
        Changed: Functions <code class="class">VectorTools</code>::<code
        class="member">point_value</code> and <code class="class">VectorTools</code>::<code

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.