OK, so I've got a simple version of what we've discussed working.
I've attached a patch, as well as some new files for the instantiation
of FEScalar.
Some notes:
- As a first cut, I used lots of hacks in dof_map.h that will only work
in serial
- I used Roy's idea of letting the Order of SCALAR indicate how many
scalar variables we want -- I guess then it would be good to make sure
that we only add one SCALAR to a system, but I haven't done that yet...
- I didn't (yet) make any changes to the sparsity pattern, so adding to
the last row(s) and column(s) of the matrix will be slow
- Plotting the SCALAR variables gives garbage at the moment, although I
guess it should be easy to change GMVIO etc to plot the constant value
of the SCALAR
I'll attach an example program to another email since if I attach it
here it puts me over the 40kb limit.
- Dave
Roy Stogner wrote:
On Sat, 18 Apr 2009, David Knezevic wrote:
So, I'd really like to add this system.add_variable(foo,SCALAR)
functionality to libMesh. I've been poking around the library, but I'm
not really sure how to get started, so I could definitely use some
pointers...
It seems to me that the important details for this are in DofMap, for
storing and retrieving the dof index of the scalar variable, as well as
for setting the sparsity pattern/number of non-zeros...?
That's right.
Some thoughts I've had are:
- Add a new enum SCALAR to FEFamily (or alternatively, a new enum SCALAR
to Order)?
FEFamily. The Order would then specify how many scalars to add.
Maybe next write an FE<SCALAR> specialization? Something that acted
sort of like a discontinuous monomial element? That might help
avoiding the need to add special cases to a lot of loops in the
library later.
- Short circuit all the loops over elements for DOF counting in DofMap
for SCALAR variables,
If the FE<SCALAR> claimed to have zero dofs per element, DofMap would
basically skip right over it... but I'm not sure that's an intuitive
behavior; this might be a loop you really do want a special case on.
and instead store the dof index of each SCALAR variable in a vector
in DofMap?
Right.
- However, I can't see where in the code one should compute the dof
index of a SCALAR variable in the first place?
That's all in dof_map.C.
- Ben, regarding setting the nonzero count and the number of rows that
you mentioned; where is this controlled? In SparsityPattern?
In SparseMatrix, I believe, although the code to build the sparsity
pattern is in DofMap.
---
Roy
Index: include/enums/enum_fe_family.h
===================================================================
--- include/enums/enum_fe_family.h (revision 3361)
+++ include/enums/enum_fe_family.h (working copy)
@@ -53,6 +53,10 @@
// C1 elements
CLOUGH = 21,
HERMITE = 22,
+
+ // A scalar variable that couples to
+ // all other DOFs in the system
+ SCALAR = 31,
INVALID_FE = 42};
Index: include/base/dof_map.h
===================================================================
--- include/base/dof_map.h (revision 3361)
+++ include/base/dof_map.h (working copy)
@@ -338,19 +338,27 @@
/**
* @returns the total number of degrees of freedom in the problem.
*/
- unsigned int n_dofs() const { return _n_dfs; }
+ unsigned int n_dofs() const { return _n_dfs + _n_SCALAR_dofs; }
/**
+ * @returns the number of SCALAR dofs.
+ */
+ unsigned int n_SCALAR_dofs() const { return _n_SCALAR_dofs; }
+
+ /**
* @returns the number of degrees of freedom on this processor.
*/
unsigned int n_local_dofs () const
- { return this->n_dofs_on_processor (libMesh::processor_id()); }
+ { return this->n_dofs_on_processor (libMesh::processor_id()) +
_n_SCALAR_dofs; }
/**
* Returns the number of degrees of freedom on subdomain \p proc.
+ *
+ * This _n_SCALAR_dofs hack won't work in parallel, just test it out.
+ * Should store the number of SCALAR dofs on each processor...
*/
unsigned int n_dofs_on_processor(const unsigned int proc) const
- { libmesh_assert(proc < _first_df.size()); return (_end_df[proc] -
_first_df[proc]); }
+ { libmesh_assert(proc < _first_df.size()); return (_end_df[proc] -
_first_df[proc] + _n_SCALAR_dofs); }
/**
* Returns the first dof index that is local to subdomain \p proc.
@@ -362,16 +370,20 @@
* Returns the last dof index that is local to subdomain \p proc.
* This function is now deprecated, because it returns nonsense in the rare
* case where \p proc has no local dof indices. Use end_dof() instead.
+ *
+ * This _n_SCALAR_dofs hack won't work in parallel, just for testing...
*/
unsigned int last_dof(const unsigned int proc = libMesh::processor_id())
const
- { libmesh_assert(proc < _end_df.size()); return (_end_df[proc] - 1); }
+ { libmesh_assert(proc < _end_df.size()); return (_end_df[proc] - 1 +
_n_SCALAR_dofs); }
/**
* Returns the first dof index that is after all indices local to subdomain
\p proc.
* Analogous to the end() member function of STL containers.
+ *
+ * This _n_SCALAR_dofs hack won't work in parallel, just for testing...
*/
unsigned int end_dof(const unsigned int proc = libMesh::processor_id()) const
- { libmesh_assert(proc < _end_df.size()); return _end_df[proc]; }
+ { libmesh_assert(proc < _end_df.size()); return _end_df[proc] +
_n_SCALAR_dofs; }
/**
* Returns the first local degree of freedom index for variable \p var.
@@ -386,12 +398,15 @@
/**
* Returns (one past) the last local degree of freedom index for variable \p
var.
* Analogous to the end() member function of STL containers.
+ *
+ * This _n_SCALAR_dofs hack won't work in parallel, just test it out.
+ * Should store the number of SCALAR dofs on each processor...
*/
unsigned int variable_last_local_dof (const unsigned int var) const
{
libmesh_assert ((var+1) < _var_first_local_df.size());
libmesh_assert (_var_first_local_df[var+1] != DofObject::invalid_id);
- return _var_first_local_df[var+1];
+ return _var_first_local_df[var+1] + _n_SCALAR_dofs;
}
/**
@@ -810,6 +825,11 @@
*/
unsigned int _n_dfs;
+ /**
+ * The number of SCALAR dofs.
+ */
+ unsigned int _n_SCALAR_dofs;
+
#ifdef LIBMESH_ENABLE_AMR
/**
@@ -850,7 +870,8 @@
_dof_coupling(NULL),
_sys_number(number),
// _matrix(NULL),
- _n_dfs(0)
+ _n_dfs(0),
+ _n_SCALAR_dofs(0)
#ifdef LIBMESH_ENABLE_AMR
, _n_old_dfs(0)
#endif
Index: include/fe/fe.h
===================================================================
--- include/fe/fe.h (revision 3361)
+++ include/fe/fe.h (working copy)
@@ -573,7 +573,23 @@
};
+//-------------------------------------------------------------
+// FEScalar class definition
+template <unsigned int Dim>
+class FEScalar : public FE<Dim,SCALAR>
+{
+public:
+ /**
+ * Constructor. Creates a SCALAR finite element
+ * which simply represents one or more
+ * extra DOFs coupled to all other DOFs in
+ * the system.
+ */
+ FEScalar(const FEType& fet);
+};
+
+
/**
* XYZ finite elements. These require specialization
* because the shape functions are defined in terms of
@@ -647,6 +663,7 @@
+
/**
* Provide Typedefs for various element types.
*/
@@ -715,6 +732,7 @@
* Monomial finite element.
*/
typedef FE<3,MONOMIAL> FEMonomial3D;
+
}
@@ -794,7 +812,7 @@
// ------------------------------------------------------------
-// FELagrange class inline members
+// FEXYZ class inline members
template <unsigned int Dim>
inline
FEXYZ<Dim>::FEXYZ (const FEType& fet) :
@@ -802,5 +820,13 @@
{
}
+// ------------------------------------------------------------
+// FEScalar class inline members
+template <unsigned int Dim>
+inline
+FEScalar<Dim>::FEScalar (const FEType& fet) :
+ FE<Dim,SCALAR> (fet)
+{
+}
#endif
Index: include/fe/fe_macro.h
===================================================================
--- include/fe/fe_macro.h (revision 3361)
+++ include/fe/fe_macro.h (working copy)
@@ -38,6 +38,7 @@
template class FE< (_dim), HIERARCHIC>; \
template class FE< (_dim), LAGRANGE>; \
template class FE< (_dim), MONOMIAL>; \
+ template class FE< (_dim), SCALAR>; \
template class FE< (_dim), XYZ>
#define INSTANTIATE_IMAP(_dim) \
@@ -46,6 +47,7 @@
template void FE<_dim,HIERARCHIC>::inverse_map(const Elem*,const
std::vector<Point>&,std::vector<Point>&,Real,bool);\
template void FE<_dim,LAGRANGE>::inverse_map(const Elem*,const
std::vector<Point>&,std::vector<Point>&,Real,bool);\
template void FE<_dim,MONOMIAL>::inverse_map(const Elem*,const
std::vector<Point>&,std::vector<Point>&,Real,bool);\
+ template void FE<_dim,SCALAR>::inverse_map(const Elem*,const
std::vector<Point>&,std::vector<Point>&,Real,bool);\
template void FE<_dim,XYZ>::inverse_map(const Elem*,const
std::vector<Point>&,std::vector<Point>&,Real,bool)
#else //LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
@@ -55,6 +57,7 @@
template class FE< (_dim), HIERARCHIC>; \
template class FE< (_dim), LAGRANGE>; \
template class FE< (_dim), MONOMIAL>; \
+ template class FE< (_dim), SCALAR>; \
template class FE< (_dim), BERNSTEIN>; \
template class FE< (_dim), SZABAB>; \
template class FE< (_dim), XYZ>
@@ -65,6 +68,7 @@
template void FE<_dim,HIERARCHIC>::inverse_map(const Elem*,const
std::vector<Point>&,std::vector<Point>&,Real,bool);\
template void FE<_dim,LAGRANGE>::inverse_map(const Elem*,const
std::vector<Point>&,std::vector<Point>&,Real,bool);\
template void FE<_dim,MONOMIAL>::inverse_map(const Elem*,const
std::vector<Point>&,std::vector<Point>&,Real,bool);\
+ template void FE<_dim,SCALAR>::inverse_map(const Elem*,const
std::vector<Point>&,std::vector<Point>&,Real,bool);\
template void FE<_dim,BERNSTEIN>::inverse_map(const Elem*,const
std::vector<Point>&,std::vector<Point>&,Real,bool);\
template void FE<_dim,SZABAB>::inverse_map(const Elem*,const
std::vector<Point>&,std::vector<Point>&,Real,bool);\
template void FE<_dim,XYZ>::inverse_map(const Elem*,const
std::vector<Point>&,std::vector<Point>&,Real,bool)
Index: src/utils/string_to_enum.C
===================================================================
--- src/utils/string_to_enum.C (revision 3361)
+++ src/utils/string_to_enum.C (working copy)
@@ -232,6 +232,7 @@
fefamily_to_enum["LAGRANGE" ]=LAGRANGE;
fefamily_to_enum["HIERARCHIC" ]=HIERARCHIC;
fefamily_to_enum["MONOMIAL" ]=MONOMIAL;
+ fefamily_to_enum["SCALAR" ]=SCALAR;
fefamily_to_enum["XYZ" ]=XYZ;
fefamily_to_enum["BERNSTEIN" ]=BERNSTEIN;
fefamily_to_enum["SZABAB" ]=SZABAB;
Index: src/base/dof_map.C
===================================================================
--- src/base/dof_map.C (revision 3361)
+++ src/base/dof_map.C (working copy)
@@ -349,15 +349,25 @@
for ( ; elem_it != elem_end; ++elem_it)
(*elem_it)->set_n_vars(this->sys_number(),n_var);
}
-
+
+ // Zero _n_SCALAR_dofs, it will be updated below.
+ this->_n_SCALAR_dofs = 0;
+
//------------------------------------------------------------
// Next allocate space for the DOF indices
for (unsigned int var=0; var<this->n_variables(); var++)
{
const System::Variable &var_description = this->variable(var);
const FEType& base_fe_type = this->variable_type(var);
-
+
+ // Skip dof counting for a SCALAR variable
+ if(base_fe_type.family == SCALAR)
+ {
+ this->_n_SCALAR_dofs = base_fe_type.order;
+ continue;
+ }
+
// This should be constant even on p-refined elements
const bool extra_hanging_dofs =
FEInterface::extra_hanging_dofs(base_fe_type);
@@ -753,6 +763,16 @@
std::fill (_var_first_local_df.begin(),
_var_first_local_df.end(),
DofObject::invalid_id);
+
+ // Set _n_SCALAR_dofs
+ this->_n_SCALAR_dofs = 0;
+ for (unsigned var=0; var<n_vars; var++)
+ {
+ if( this->variable(var).type().family != SCALAR )
+ {
+ _n_SCALAR_dofs = this->variable(var).type().order;
+ }
+ }
//-------------------------------------------------------------------------
// First count and assign temporary numbers to local dofs
@@ -770,9 +790,11 @@
for (unsigned int n=0; n<n_nodes; n++)
{
Node* node = elem->get_node(n);
-
+
for (unsigned var=0; var<n_vars; var++)
- if (this->variable(var).active_on_subdomain(elem->subdomain_id()))
+ {
+ if( (this->variable(var).type().family != SCALAR) &&
+
(this->variable(var).active_on_subdomain(elem->subdomain_id())) )
{
// assign dof numbers (all at once) if this is
// our node and if they aren't already there
@@ -788,11 +810,13 @@
next_free_dof += node->n_comp(sys_num,var);
}
}
+ }
}
-
+
// Now number the element DOFS
for (unsigned var=0; var<n_vars; var++)
- if (this->variable(var).active_on_subdomain(elem->subdomain_id()))
+ if ( (this->variable(var).type().family != SCALAR) &&
+ (this->variable(var).active_on_subdomain(elem->subdomain_id())) )
if (elem->n_comp(sys_num,var) > 0)
{
libmesh_assert (elem->dof_number(sys_num,var,0) ==
@@ -843,15 +867,25 @@
// We will cache the first local index for each variable
_var_first_local_df.clear();
+
+ // Zero _n_SCALAR_dofs, we will set it below
+ this->_n_SCALAR_dofs = 0;
//-------------------------------------------------------------------------
// First count and assign temporary numbers to local dofs
for (unsigned var=0; var<n_vars; var++)
{
_var_first_local_df.push_back(next_free_dof);
-
+
const System::Variable var_description = this->variable(var);
+ // Skip dof counting for a SCALAR variable
+ if(var_description.type().family == SCALAR)
+ {
+ this->_n_SCALAR_dofs = var_description.type().order;
+ continue;
+ }
+
MeshBase::element_iterator elem_it =
mesh.active_local_elements_begin();
const MeshBase::element_iterator elem_end =
mesh.active_local_elements_end();
@@ -1249,9 +1283,24 @@
unsigned int tot_size = 0;
#endif
+ bool SCALAR_var_flag = false;
+
// Get the dof numbers
for (unsigned int v=0; v<n_vars; v++)
if ((v == vn) || (vn == libMesh::invalid_uint))
+ {
+ if(this->variable(v).type().family == SCALAR)
+ {
+ SCALAR_var_flag = true;
+
+#ifdef DEBUG
+ FEType fe_type = this->variable_type(v);
+ tot_size += FEInterface::n_dofs(dim,
+ fe_type,
+ type);
+#endif
+ }
+ else
if (this->variable(v).active_on_subdomain(elem->subdomain_id()))
{ // Do this for all the variables if one was not specified
// or just for the specified variable
@@ -1346,8 +1395,15 @@
di.resize(di.size() + nc, DofObject::invalid_id);
}
}
- } // end loop over variables
+ }
+ } // end loop over variables
+ // Finally append any SCALAR dofs if we asked for them
+ // or no variable number was specified
+ if ( SCALAR_var_flag || (vn == libMesh::invalid_uint) )
+ for(unsigned int i = (n_dofs() - n_SCALAR_dofs()); i < n_dofs(); i++)
+ di.push_back(i);
+
#ifdef DEBUG
libmesh_assert (tot_size == di.size());
#endif
@@ -1386,7 +1442,8 @@
// Get the dof numbers
for (unsigned int v=0; v<n_vars; v++)
if ((v == vn) || (vn == libMesh::invalid_uint))
- if (this->variable(v).active_on_subdomain(elem->subdomain_id()))
+ if ( (this->variable(v).type().family != SCALAR) &&
+ (this->variable(v).active_on_subdomain(elem->subdomain_id())) )
{ // Do this for all the variables if one was not specified
// or just for the specified variable
Index: src/fe/fe_interface.C
===================================================================
--- src/fe/fe_interface.C (revision 3361)
+++ src/fe/fe_interface.C (working copy)
@@ -79,6 +79,9 @@
case MONOMIAL:
return FE<1,MONOMIAL>::n_shape_functions(t, o);
+ case SCALAR:
+ return FE<1,SCALAR>::n_shape_functions(t, o);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -119,6 +122,9 @@
case MONOMIAL:
return FE<2,MONOMIAL>::n_shape_functions(t, o);
+ case SCALAR:
+ return FE<2,SCALAR>::n_shape_functions(t, o);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -155,6 +161,9 @@
case MONOMIAL:
return FE<3,MONOMIAL>::n_shape_functions(t, o);
+ case SCALAR:
+ return FE<3,SCALAR>::n_shape_functions(t, o);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -222,6 +231,9 @@
case MONOMIAL:
return FE<1,MONOMIAL>::n_dofs(t, o);
+ case SCALAR:
+ return FE<1,SCALAR>::n_dofs(t, o);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -261,6 +273,9 @@
case MONOMIAL:
return FE<2,MONOMIAL>::n_dofs(t, o);
+ case SCALAR:
+ return FE<2,SCALAR>::n_dofs(t, o);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -297,6 +312,9 @@
case MONOMIAL:
return FE<3,MONOMIAL>::n_dofs(t, o);
+ case SCALAR:
+ return FE<3,SCALAR>::n_dofs(t, o);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -364,6 +382,9 @@
case MONOMIAL:
return FE<1,MONOMIAL>::n_dofs_at_node(t, o, n);
+ case SCALAR:
+ return FE<1,SCALAR>::n_dofs_at_node(t, o, n);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -403,6 +424,9 @@
case MONOMIAL:
return FE<2,MONOMIAL>::n_dofs_at_node(t, o, n);
+ case SCALAR:
+ return FE<2,SCALAR>::n_dofs_at_node(t, o, n);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -439,6 +463,9 @@
case MONOMIAL:
return FE<3,MONOMIAL>::n_dofs_at_node(t, o, n);
+ case SCALAR:
+ return FE<3,SCALAR>::n_dofs_at_node(t, o, n);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -506,6 +533,9 @@
case MONOMIAL:
return FE<1,MONOMIAL>::n_dofs_per_elem(t, o);
+ case SCALAR:
+ return FE<1,SCALAR>::n_dofs_per_elem(t, o);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -545,6 +575,9 @@
case MONOMIAL:
return FE<2,MONOMIAL>::n_dofs_per_elem(t, o);
+ case SCALAR:
+ return FE<2,SCALAR>::n_dofs_per_elem(t, o);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -581,6 +614,9 @@
case MONOMIAL:
return FE<3,MONOMIAL>::n_dofs_per_elem(t, o);
+ case SCALAR:
+ return FE<3,SCALAR>::n_dofs_per_elem(t, o);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -647,6 +683,10 @@
FE<1,MONOMIAL>::dofs_on_side(elem, o, s, di);
return;
+ case SCALAR:
+ FE<1,SCALAR>::dofs_on_side(elem, o, s, di);
+ return;
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -694,6 +734,10 @@
FE<2,MONOMIAL>::dofs_on_side(elem, o, s, di);
return;
+ case SCALAR:
+ FE<2,SCALAR>::dofs_on_side(elem, o, s, di);
+ return;
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -737,6 +781,10 @@
FE<3,MONOMIAL>::dofs_on_side(elem, o, s, di);
return;
+ case SCALAR:
+ FE<3,SCALAR>::dofs_on_side(elem, o, s, di);
+ return;
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -803,6 +851,10 @@
FE<1,MONOMIAL>::dofs_on_edge(elem, o, e, di);
return;
+ case SCALAR:
+ FE<1,SCALAR>::dofs_on_edge(elem, o, e, di);
+ return;
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -850,6 +902,10 @@
FE<2,MONOMIAL>::dofs_on_edge(elem, o, e, di);
return;
+ case SCALAR:
+ FE<2,SCALAR>::dofs_on_edge(elem, o, e, di);
+ return;
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -893,6 +949,10 @@
FE<3,MONOMIAL>::dofs_on_edge(elem, o, e, di);
return;
+ case SCALAR:
+ FE<3,SCALAR>::dofs_on_edge(elem, o, e, di);
+ return;
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -975,6 +1035,11 @@
elem_soln, nodal_soln);
return;
+ case SCALAR:
+ FE<1,SCALAR>::nodal_soln(elem, order,
+ elem_soln, nodal_soln);
+ return;
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -1030,6 +1095,11 @@
elem_soln, nodal_soln);
return;
+ case SCALAR:
+ FE<2,SCALAR>::nodal_soln(elem, order,
+ elem_soln, nodal_soln);
+ return;
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -1080,6 +1150,11 @@
elem_soln, nodal_soln);
return;
+ case SCALAR:
+ FE<3,SCALAR>::nodal_soln(elem, order,
+ elem_soln, nodal_soln);
+ return;
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -1149,6 +1224,9 @@
case MONOMIAL:
return FE<1,MONOMIAL>::inverse_map(elem, p, tolerance, secure);
+ case SCALAR:
+ return FE<1,SCALAR>::inverse_map(elem, p, tolerance, secure);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -1188,6 +1266,9 @@
case MONOMIAL:
return FE<2,MONOMIAL>::inverse_map(elem, p, tolerance, secure);
+ case SCALAR:
+ return FE<2,SCALAR>::inverse_map(elem, p, tolerance, secure);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -1223,6 +1304,9 @@
case MONOMIAL:
return FE<3,MONOMIAL>::inverse_map(elem, p, tolerance, secure);
+ case SCALAR:
+ return FE<3,SCALAR>::inverse_map(elem, p, tolerance, secure);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -1319,6 +1403,10 @@
case MONOMIAL:
FE<1,MONOMIAL>::inverse_map(elem, physical_points,
reference_points, tolerance, secure);
return;
+
+ case SCALAR:
+ FE<1,SCALAR>::inverse_map(elem, physical_points, reference_points,
tolerance, secure);
+ return;
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
@@ -1366,6 +1454,10 @@
case MONOMIAL:
FE<2,MONOMIAL>::inverse_map(elem, physical_points,
reference_points, tolerance, secure);
return;
+
+ case SCALAR:
+ FE<2,SCALAR>::inverse_map(elem, physical_points, reference_points,
tolerance, secure);
+ return;
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
@@ -1409,6 +1501,10 @@
case MONOMIAL:
FE<3,MONOMIAL>::inverse_map(elem, physical_points,
reference_points, tolerance, secure);
return;
+
+ case SCALAR:
+ FE<3,SCALAR>::inverse_map(elem, physical_points, reference_points,
tolerance, secure);
+ return;
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
@@ -1490,6 +1586,9 @@
case MONOMIAL:
return FE<1,MONOMIAL>::shape(t,o,i,p);
+ case SCALAR:
+ return FE<1,SCALAR>::shape(t,o,i,p);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -1529,6 +1628,9 @@
case MONOMIAL:
return FE<2,MONOMIAL>::shape(t,o,i,p);
+ case SCALAR:
+ return FE<2,SCALAR>::shape(t,o,i,p);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -1565,6 +1667,9 @@
case MONOMIAL:
return FE<3,MONOMIAL>::shape(t,o,i,p);
+ case SCALAR:
+ return FE<3,SCALAR>::shape(t,o,i,p);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -1633,6 +1738,9 @@
case MONOMIAL:
return FE<1,MONOMIAL>::shape(elem,o,i,p);
+ case SCALAR:
+ return FE<1,SCALAR>::shape(elem,o,i,p);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -1672,6 +1780,9 @@
case MONOMIAL:
return FE<2,MONOMIAL>::shape(elem,o,i,p);
+ case SCALAR:
+ return FE<2,SCALAR>::shape(elem,o,i,p);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
@@ -1708,6 +1819,9 @@
case MONOMIAL:
return FE<3,MONOMIAL>::shape(elem,o,i,p);
+ case SCALAR:
+ return FE<3,SCALAR>::shape(elem,o,i,p);
+
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
case BERNSTEIN:
Index: src/fe/fe_base.C
===================================================================
--- src/fe/fe_base.C (revision 3361)
+++ src/fe/fe_base.C (working copy)
@@ -103,6 +103,12 @@
return ap;
}
+ case SCALAR:
+ {
+ AutoPtr<FEBase> ap(new FEScalar<0>(fet));
+ return ap;
+ }
+
default:
std::cout << "ERROR: Bad FEType.family= " << fet.family <<
std::endl;
libmesh_error();
@@ -163,6 +169,12 @@
return ap;
}
+ case SCALAR:
+ {
+ AutoPtr<FEBase> ap(new FEScalar<1>(fet));
+ return ap;
+ }
+
default:
std::cout << "ERROR: Bad FEType.family= " << fet.family <<
std::endl;
libmesh_error();
@@ -225,6 +237,12 @@
return ap;
}
+ case SCALAR:
+ {
+ AutoPtr<FEBase> ap(new FEScalar<2>(fet));
+ return ap;
+ }
+
default:
std::cout << "ERROR: Bad FEType.family= " << fet.family <<
std::endl;
libmesh_error();
@@ -288,6 +306,12 @@
return ap;
}
+ case SCALAR:
+ {
+ AutoPtr<FEBase> ap(new FEScalar<3>(fet));
+ return ap;
+ }
+
default:
std::cout << "ERROR: Bad FEType.family= " << fet.family <<
std::endl;
libmesh_error();
Index: src/fe/fe_boundary.C
===================================================================
--- src/fe/fe_boundary.C (revision 3361)
+++ src/fe/fe_boundary.C (working copy)
@@ -54,6 +54,8 @@
REINIT_ERROR(0, LAGRANGE, edge_reinit)
REINIT_ERROR(0, MONOMIAL, reinit)
REINIT_ERROR(0, MONOMIAL, edge_reinit)
+REINIT_ERROR(0, SCALAR, reinit)
+REINIT_ERROR(0, SCALAR, edge_reinit)
REINIT_ERROR(0, XYZ, reinit)
REINIT_ERROR(0, XYZ, edge_reinit)
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
@@ -69,6 +71,7 @@
REINIT_ERROR(1, LAGRANGE, edge_reinit)
REINIT_ERROR(1, XYZ, edge_reinit)
REINIT_ERROR(1, MONOMIAL, edge_reinit)
+REINIT_ERROR(1, SCALAR, edge_reinit)
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
REINIT_ERROR(1, BERNSTEIN, edge_reinit)
REINIT_ERROR(1, SZABAB, edge_reinit)
@@ -706,6 +709,7 @@
template void FE<1,CLOUGH>::reinit(Elem const*, unsigned int, Real);
template void FE<1,HERMITE>::reinit(Elem const*, unsigned int, Real);
template void FE<1,MONOMIAL>::reinit(Elem const*, unsigned int, Real);
+template void FE<1,SCALAR>::reinit(Elem const*, unsigned int, Real);
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
template void FE<1,BERNSTEIN>::reinit(Elem const*, unsigned int, Real);
template void FE<1,SZABAB>::reinit(Elem const*, unsigned int, Real);
@@ -722,6 +726,8 @@
template void FE<2,HERMITE>::edge_reinit(Elem const*, unsigned int, Real);
template void FE<2,MONOMIAL>::reinit(Elem const*, unsigned int, Real);
template void FE<2,MONOMIAL>::edge_reinit(Elem const*, unsigned int, Real);
+template void FE<2,SCALAR>::reinit(Elem const*, unsigned int, Real);
+template void FE<2,SCALAR>::edge_reinit(Elem const*, unsigned int, Real);
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
template void FE<2,BERNSTEIN>::reinit(Elem const*, unsigned int, Real);
template void FE<2,BERNSTEIN>::edge_reinit(Elem const*, unsigned int, Real);
@@ -744,6 +750,8 @@
template void FE<3,HERMITE>::edge_reinit(Elem const*, unsigned int, Real);
template void FE<3,MONOMIAL>::reinit(Elem const*, unsigned int, Real);
template void FE<3,MONOMIAL>::edge_reinit(Elem const*, unsigned int, Real);
+template void FE<3,SCALAR>::reinit(Elem const*, unsigned int, Real);
+template void FE<3,SCALAR>::edge_reinit(Elem const*, unsigned int, Real);
#ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
template void FE<3,BERNSTEIN>::reinit(Elem const*, unsigned int, Real);
template void FE<3,BERNSTEIN>::edge_reinit(Elem const*, unsigned int, Real);
// $Id: fe_lagrange.C 3244 2009-01-30 23:42:03Z roystgnr $
// The libMesh Finite Element Library.
// Copyright (C) 2002-2008 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
// Lesser General Public License for more details.
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
// Local includes
#include "dof_map.h"
#include "fe.h"
#include "fe_macro.h"
#include "elem.h"
// ------------------------------------------------------------
// SCALAR-specific implementations
template <unsigned int Dim, FEFamily T>
void FE<Dim,T>::nodal_soln(const Elem*,
const Order,
const std::vector<Number>&,
std::vector<Number>&)
{
// Do nothing, not relevant to SCALAR
return;
}
template <unsigned int Dim, FEFamily T>
unsigned int FE<Dim,T>::n_dofs(const ElemType, const Order o)
{
// The Order indicates the number of SCALAR dofs
return o;
}
template <unsigned int Dim, FEFamily T>
unsigned int FE<Dim,T>::n_dofs_at_node(const ElemType,
const Order,
const unsigned int)
{
// SCALARs have no dofs at nodes
return 0;
}
template <unsigned int Dim, FEFamily T>
unsigned int FE<Dim,T>::n_dofs_per_elem(const ElemType,
const Order)
{
// SCALARs have no dofs per element
return 0;
}
template <unsigned int Dim, FEFamily T>
FEContinuity FE<Dim,T>::get_continuity() const
{
// This doesn't really make sense for a SCALAR...
return C_ZERO;
}
template <unsigned int Dim, FEFamily T>
bool FE<Dim,T>::is_hierarchic() const
{
return false;
}
#ifdef LIBMESH_ENABLE_AMR
template <unsigned int Dim, FEFamily T>
void FE<Dim,T>::compute_constraints (DofConstraints &constraints,
DofMap &dof_map,
const unsigned int variable_number,
const Elem* elem)
{
return;
}
#endif // #ifdef LIBMESH_ENABLE_AMR
template <unsigned int Dim, FEFamily T>
bool FE<Dim,T>::shapes_need_reinit() const
{
return false;
}
//--------------------------------------------------------------
// Explicit instantiation of member functions
INSTANTIATE_MBRF(0,SCALAR);
INSTANTIATE_MBRF(1,SCALAR);
INSTANTIATE_MBRF(2,SCALAR);
INSTANTIATE_MBRF(3,SCALAR);// $Id: fe_lagrange_shape_0D.C 3243 2009-01-30 23:04:41Z roystgnr $
// The libMesh Finite Element Library.
// Copyright (C) 2002-2008 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
// Lesser General Public License for more details.
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
// C++ includes
// Local includes
#include "fe.h"
#include "elem.h"
template <>
Real FE<0,SCALAR>::shape(const ElemType,
const Order,
const unsigned int,
const Point&)
{
libmesh_error();
return 1.;
}
template <>
Real FE<0,SCALAR>::shape(const Elem*,
const Order,
const unsigned int,
const Point&)
{
libmesh_error();
return 1.;
}
template <>
Real FE<0,SCALAR>::shape_deriv(const ElemType,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}
template <>
Real FE<0,SCALAR>::shape_deriv(const Elem*,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}
template <>
Real FE<0,SCALAR>::shape_second_deriv(const ElemType,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}
template <>
Real FE<0,SCALAR>::shape_second_deriv(const Elem*,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}// $Id: fe_lagrange_shape_0D.C 3243 2009-01-30 23:04:41Z roystgnr $
// The libMesh Finite Element Library.
// Copyright (C) 2002-2008 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
// Lesser General Public License for more details.
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
// C++ includes
// Local includes
#include "fe.h"
#include "elem.h"
template <>
Real FE<1,SCALAR>::shape(const ElemType,
const Order,
const unsigned int,
const Point&)
{
libmesh_error();
return 1.;
}
template <>
Real FE<1,SCALAR>::shape(const Elem*,
const Order,
const unsigned int,
const Point&)
{
libmesh_error();
return 1.;
}
template <>
Real FE<1,SCALAR>::shape_deriv(const ElemType,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}
template <>
Real FE<1,SCALAR>::shape_deriv(const Elem*,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}
template <>
Real FE<1,SCALAR>::shape_second_deriv(const ElemType,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}
template <>
Real FE<1,SCALAR>::shape_second_deriv(const Elem*,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}// $Id: fe_lagrange_shape_0D.C 3243 2009-01-30 23:04:41Z roystgnr $
// The libMesh Finite Element Library.
// Copyright (C) 2002-2008 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
// Lesser General Public License for more details.
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
// C++ includes
// Local includes
#include "fe.h"
#include "elem.h"
template <>
Real FE<2,SCALAR>::shape(const ElemType,
const Order,
const unsigned int,
const Point&)
{
libmesh_error();
return 1.;
}
template <>
Real FE<2,SCALAR>::shape(const Elem*,
const Order,
const unsigned int,
const Point&)
{
libmesh_error();
return 1.;
}
template <>
Real FE<2,SCALAR>::shape_deriv(const ElemType,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}
template <>
Real FE<2,SCALAR>::shape_deriv(const Elem*,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}
template <>
Real FE<2,SCALAR>::shape_second_deriv(const ElemType,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}
template <>
Real FE<2,SCALAR>::shape_second_deriv(const Elem*,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}// $Id: fe_lagrange_shape_0D.C 3243 2009-01-30 23:04:41Z roystgnr $
// The libMesh Finite Element Library.
// Copyright (C) 2002-2008 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
// Lesser General Public License for more details.
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
// C++ includes
// Local includes
#include "fe.h"
#include "elem.h"
template <>
Real FE<3,SCALAR>::shape(const ElemType,
const Order,
const unsigned int,
const Point&)
{
libmesh_error();
return 1.;
}
template <>
Real FE<3,SCALAR>::shape(const Elem*,
const Order,
const unsigned int,
const Point&)
{
libmesh_error();
return 1.;
}
template <>
Real FE<3,SCALAR>::shape_deriv(const ElemType,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}
template <>
Real FE<3,SCALAR>::shape_deriv(const Elem*,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}
template <>
Real FE<3,SCALAR>::shape_second_deriv(const ElemType,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}
template <>
Real FE<3,SCALAR>::shape_second_deriv(const Elem*,
const Order,
const unsigned int,
const unsigned int,
const Point&)
{
libmesh_error();
return 0.;
}------------------------------------------------------------------------------
The NEW KODAK i700 Series Scanners deliver under ANY circumstances! Your
production scanning environment may not be a perfect world - but thanks to
Kodak, there's a perfect scanner to get the job done! With the NEW KODAK i700
Series Scanner you'll get full speed at 300 dpi even with all image
processing features enabled. http://p.sf.net/sfu/kodak-com
_______________________________________________
Libmesh-devel mailing list
Libmesh-devel@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/libmesh-devel
