From: Wolfgang Bangerth Date: Mon, 18 May 1998 08:40:34 +0000 (+0000) Subject: Implement interpolation. X-Git-Tag: v8.0.0~22950 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=e43eb4101b939dcfc0801023492fa1017fc0ee8d;p=dealii.git Implement interpolation. git-svn-id: https://svn.dealii.org/trunk@304 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/numerics/vectors.cc b/deal.II/deal.II/source/numerics/vectors.cc index 040410d3c7..061aa87c06 100644 --- a/deal.II/deal.II/source/numerics/vectors.cc +++ b/deal.II/deal.II/source/numerics/vectors.cc @@ -3,6 +3,8 @@ #include #include +#include +#include #include #include #include @@ -18,9 +20,29 @@ template void VectorCreator::interpolate (const DoFHandler &dof, const FiniteElement &fe, + const Boundary &boundary, const Function &function, dVector &vec) { - ; + DoFHandler::active_cell_iterator cell = dof.begin_active(), + endc = dof.end(); + vector dofs_on_cell (fe.total_dofs); + vector dof_values_on_cell (fe.total_dofs); + vector > ansatz_points (fe.total_dofs); + for (; cell!=endc; ++cell) + { + // for each cell: + // get location of finite element + // off-points + fe.get_ansatz_points (cell, boundary, ansatz_points); + // get function values at these points + function.value_list (ansatz_points, dof_values_on_cell); + // get indices of the dofs on this cell + cell->get_dof_indices (dofs_on_cell); + // distribute function values to the + // whole vector + for (unsigned int i=0; i::project (const DoFHandler &dof, dof.n_dofs(), dof.max_couplings_between_dofs()); dof.make_sparsity_pattern (sparsity); - + constraints.condense (sparsity); + dSMatrix mass_matrix (sparsity); dVector tmp (mass_matrix.n()); MatrixCreator::create_mass_matrix (dof, fe, q, boundary, - mass_matrix, function, vec); + mass_matrix, function, tmp); + + constraints.condense (mass_matrix); + constraints.condense (tmp); int max_iter = 4000; double tolerance = 1.e-16;