Skip to content
Snippets Groups Projects
Commit 9b35bf98 authored by Oliver Sander's avatar Oliver Sander Committed by sander
Browse files

Make CubeCharacteristicFunctional inherit from NonsmoothConvexFunctional

[[Imported from SVN: r12811]]
parent 7826147e
Branches
No related tags found
No related merge requests found
#ifndef DUNE_TNNMG_FUNCTIONALS_CUBECHARACTERISTICFUNCTIONAL_HH #ifndef DUNE_TNNMG_FUNCTIONALS_CUBECHARACTERISTICFUNCTIONAL_HH
#define DUNE_TNNMG_FUNCTIONALS_CUBECHARACTERISTICFUNCTIONAL_HH #define DUNE_TNNMG_FUNCTIONALS_CUBECHARACTERISTICFUNCTIONAL_HH
#include <dune/tnnmg/problem-classes/pointwisenonlinearity.hh> #include <dune/tnnmg/functionals/nonsmoothconvexfunctional.hh>
namespace Dune { namespace Dune {
...@@ -14,10 +14,10 @@ namespace Dune { ...@@ -14,10 +14,10 @@ namespace Dune {
* \tparam LocalMatrixType Matrix type for matrices R^{N \times N} * \tparam LocalMatrixType Matrix type for matrices R^{N \times N}
*/ */
template < class LocalVectorType=Dune::FieldVector<double,1>, class LocalMatrixType=Dune::FieldMatrix<double,1,1> > template < class LocalVectorType=Dune::FieldVector<double,1>, class LocalMatrixType=Dune::FieldMatrix<double,1,1> >
class CubeCharacteristicFunctional: public PointwiseNonlinearity<LocalVectorType, LocalMatrixType> class CubeCharacteristicFunctional: public NonsmoothConvexFunctional<LocalVectorType, LocalMatrixType>
{ {
public: public:
typedef Nonlinearity<LocalVectorType, LocalMatrixType> NonlinearityType; typedef NonsmoothConvexFunctional<LocalVectorType, LocalMatrixType> NonlinearityType;
using NonlinearityType::block_size; using NonlinearityType::block_size;
...@@ -92,6 +92,26 @@ namespace Dune { ...@@ -92,6 +92,26 @@ namespace Dune {
return; return;
} }
/** \brief Set the internal position vector u_pos to v.
*
* This method is empty, because it is only needed if the nonlinearity
* does not decouple in the Euclidean directions.
*/
virtual void setVector(const VectorType& v) {};
/** \brief Update the (i,j)-th entry of the internal position vector u_pos to x.
*
* This method is empty, because it is only needed if the nonlinearity
* does not decouple in the Euclidean directions.
*
* \param i global index
* \param j local index
* \param x new value of the entry (u_pos)_(i,j)
*/
virtual void updateEntry(int i, double x, int j) {};
private: private:
const double TOL; const double TOL;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment