DirichletBC.inc.hpp 1.53 KB
Newer Older
1 2 3 4
#pragma once

#include <dune/functions/functionspacebases/interpolate.hh>

5 6
#include "Log.hpp"

7 8
namespace AMDiS
{
9 10
  template <class WorldVector>
    template <class RowFeSpace, class ColFeSpace, class DOFMatrix, class DOFVector>
11 12 13 14 15 16 17 18 19
  void DirichletBC<WorldVector>::init(bool apply, 
				      RowFeSpace const& rowFeSpace, 
				      ColFeSpace const& colFeSpace, 
				      DOFMatrix* matrix, 
				      DOFVector* rhs, 
				      DOFVector* solution)
  {
    using Dune::Functions::interpolate;
    
20
    if (!initialized) {
21
      interpolate(rowFeSpace, dirichletNodes, predicate);
22 23
      initialized = true;
    }
24 25 26
  }
    
    
27 28
  template <class WorldVector>
    template <class RowFeSpace, class ColFeSpace, class DOFMatrix, class DOFVector>
29 30 31 32 33 34 35 36 37
  void DirichletBC<WorldVector>::finish(bool apply, 
					RowFeSpace const& rowFeSpace, 
					ColFeSpace const& colFeSpace, 
					DOFMatrix* matrix, 
					DOFVector* rhs, 
					DOFVector* solution)
  {
    using Dune::Functions::interpolate;
    
38
    AMDIS_TEST_EXIT( initialized, "Boundary condition not initialized!" );
39 40
    
    // loop over the matrix rows
41
    for (size_t i = 0; i < matrix->N(); ++i) {
42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57
	if (dirichletNodes[i]) {
	    auto cIt = (*matrix)[i].begin();
	    auto cEndIt = (*matrix)[i].end();
	    // loop over nonzero matrix entries in current row
	    for (; cIt != cEndIt; ++cIt)
		*cIt = (apply && i == cIt.index()) ? 1.0 : 0.0;
	}
    }

    if (apply) {
	interpolate(rowFeSpace, *rhs, values, dirichletNodes);
	interpolate(colFeSpace, *solution, values, dirichletNodes);
    }
  }

} // end namespace AMDiS