ExcMatrixNotSquare());
-// store for each line of the matrix
+ // store for each line of the matrix
// its new line number
// after compression. If the shift is
// -1, this line will be condensed away
new_line.push_back (row-shift);
-next_constraint = lines.begin();
+ next_constraint = lines.begin();
// note: in this loop we need not check
// whether #next_constraint# is a valid
// iterator, since #next_constraint# is
};
-// dofs of quad
+ // dofs of quad
if (selected_fe->dofs_per_quad > 0)
for (unsigned int d=0; d<selected_fe->dofs_per_quad; ++d)
cell->set_dof_index (d, next_free_dof++);
-// note that this cell has been processed
+ // note that this cell has been processed
cell->set_user_flag ();
return next_free_dof;
};
-// dofs of hex
+ // dofs of hex
if (selected_fe->dofs_per_hex > 0)
for (unsigned int d=0; d<selected_fe->dofs_per_hex; ++d)
cell->set_dof_index (d, next_free_dof++);
-// note that this cell has been processed
+ // note that this cell has been processed
cell->set_user_flag ();
return next_free_dof;
27*selected_fe->dofs_per_hex);
-Assert (false, ExcNotImplemented());
+ Assert (false, ExcNotImplemented());
return 0;
};
};
-// store next free dof index
+ // store next free dof index
unsigned int next_free_number = 0;
// enumerate the first round dofs
};
-// store next free dof index
+ // store next free dof index
unsigned int next_free_number = 0;
// enumerate the first round dofs
continue;
-// store for each coordination
+ // store for each coordination
// number the dofs with these
// coordination number
multimap<unsigned int, int> dofs_by_coordination;
[dof.get_fe().system_to_component_index(j).first];
-vector<unsigned int> dofs_on_this_cell(dofs_per_cell);
+ vector<unsigned int> dofs_on_this_cell(dofs_per_cell);
DoFHandler<dim>::active_cell_iterator cell = dof.begin_active(),
endc = dof.end();
for (; cell!=endc; ++cell)
cell->face(face)->set_user_flag();
-// loop over all lines; only on lines
+ // loop over all lines; only on lines
// there can be constraints.
for (line = dof_handler.begin_line(); line != endl; ++line)
// if dofs on this line are subject
cell->face(face)->set_user_flag();
-// loop over all faces; only on faces
+ // loop over all faces; only on faces
// there can be constraints.
for (face=dof_handler.begin_face(); face != endf; ++face)
// if dofs on this line are subject
ExcWrongFieldDimension(support_points.size(), dofs_per_cell));
-// local mesh width
+ // local mesh width
const double h=(cell->vertex(1)(0) - cell->vertex(0)(0));
for (unsigned int i=0; i<q_points.size(); ++i)
};
+
template <>
void FiniteElement<1>::fill_fe_face_values (const DoFHandler<1>::cell_iterator &,
const unsigned int ,
od:
-# these are the basis functions differentiated with respect to
+ # these are the basis functions differentiated with respect to
# xi and eta. we need them for the computation of the jacobi
# matrix, since we can't just differentiate a function.
phi_xi[0] := proc(x,y) if(y<1-x) then -1; else 0; fi; end:
end:
-# Now for the mass matrix: we divide the entire cell into four
+ # Now for the mass matrix: we divide the entire cell into four
# sectors:
#
# *-------*
ExcWrongFieldDimension(support_points.size(), dofs_per_cell));
-unsigned int n_points=unit_points.size();
+ unsigned int n_points=unit_points.size();
// we need the support points in any
// way, wanted or not by the user
};
-if (compute_jacobians_grad)
+ if (compute_jacobians_grad)
switch (dim)
{
case 1:
prolongation[7](63,62) = 225.0/256.0;
-restriction[0](0,0) = 1.0;
+ restriction[0](0,0) = 1.0;
restriction[0](8,9) = 1.0;
restriction[0](14,15) = 1.0;
restriction[0](24,25) = 1.0;
od;
-# to get the restriction (interpolation) matrices, evaluate
+ # to get the restriction (interpolation) matrices, evaluate
# the basis functions on the child cells at the global
# interpolation points
};
-// now set up the interface constraints.
+ // now set up the interface constraints.
// this is kind'o hairy, so don't try
// to do it dimension independent
build_interface_constraints ();
Assert((i<dofs_per_cell), ExcIndexRange(i, 0, dofs_per_cell));
-pair<unsigned,unsigned> comp = system_to_component_index(i);
+ pair<unsigned,unsigned> comp = system_to_component_index(i);
return base_element(component_to_base_table[comp.first])
.shape_grad_grad(comp.second, p);
contract (tmp2, tmp1, 1, jacobi_matrices[j], 1);
-// second part:
+ // second part:
// tmp1 := (d_k J_lj) (d_l phi)
contract (tmp1, jacobi_matrices_grad[j], 2, unit_shape_gradients[i][j]);
// tmp1_kj J_ki
};
-// compute Jacobi determinants in
+ // compute Jacobi determinants in
// quadrature points.
// refer to the general doc for
// why we take the inverse of the
};
-Tensor<2,dim> tmp1, tmp2;
+ Tensor<2,dim> tmp1, tmp2;
if (update_flags & update_second_derivatives)
for (unsigned int i=0; i<fe->dofs_per_cell; ++i)
for (unsigned int j=0; j<n_quadrature_points; ++j)
contract (tmp2, tmp1, 1, jacobi_matrices[j], 1);
-// second part:
+ // second part:
// tmp1 := (d_k J_lj) (d_l phi)
contract (tmp1,
jacobi_matrices_grad[j], 2,
};
-// compute Jacobi determinants in
+ // compute Jacobi determinants in
// quadrature points.
// refer to the general doc for
// why we take the inverse of the
contract (tmp2, tmp1, 1, jacobi_matrices[j], 1);
-// second part:
+ // second part:
// tmp1 := (d_k J_lj) (d_l phi)
contract (tmp1,
jacobi_matrices_grad[j], 2,
};
-// compute Jacobi determinants in
+ // compute Jacobi determinants in
// quadrature points.
// refer to the general doc for
// why we take the inverse of the
ExcWrongFieldDimension(support_points.size(), dofs_per_cell));
-unsigned int n_points=unit_points.size();
+ unsigned int n_points=unit_points.size();
Point<dim> vertices[GeometryInfo<dim>::vertices_per_cell];
for (unsigned int l=0; l<GeometryInfo<dim>::vertices_per_cell; ++l)
vertices[l] = cell->vertex(l);
-if (compute_q_points)
+ if (compute_q_points)
{
// initialize points to zero
for (unsigned int i=0; i<n_points; ++i)
};
-One last note regarding whether we have to invert M or M transposed: it is
+ One last note regarding whether we have to invert M or M transposed: it is
easy to try out, by computing the gradients of a function on a distorted
cell (just move one vertex) where the nodal values for linear elements
are one for the moved vertex and zero otherwise. Please also note that
fevalues.reinit (dof.begin_active(),b);
-Vector<double> val(4);
+ Vector<double> val(4);
val(2) = 1;
vector<Point<2> > grads(4);
};
-// write preamble
+ // write preamble
if (true)
{
/*
};
-DoFHandler<dim>::active_cell_iterator cell;
+ DoFHandler<dim>::active_cell_iterator cell;
DoFHandler<dim>::active_cell_iterator endc = dofs->end();
QTrapez<1> q_trapez;
};
-break;
+ break;
default:
Assert (false, ExcNotImplemented());
};
-DoFHandler<dim>::active_cell_iterator cell;
+ DoFHandler<dim>::active_cell_iterator cell;
DoFHandler<dim>::active_cell_iterator endc = dofs->end();
unsigned int cell_index=0;
out << endl << endl;
-break;
+ break;
default:
Assert (false, ExcNotImplemented());
<< endl;
-// write frame object
+ // write frame object
out << "mesh {" << endl;
DoFHandler<2>::active_cell_iterator cell;
};
-// Make output values local by
+ // Make output values local by
// copying them to a multiset.
// Perform the necessary turn.
const DoFHandler<2>::active_cell_iterator endc = dofs->end();
cells2.clear();
-// If we want shaded output we can
+ // If we want shaded output we can
// do the shading now.
if (cell_shade_p)
{
double light_norm, normal_norm;
float color;
- for (multiset<DataOut_Old<2>::EpsCellData>::iterator c=cells.begin();c!=cells.end();++c)
+ for (multiset<DataOut_Old<2>::EpsCellData>::iterator c=cells.begin();
+ c!=cells.end(); ++c)
{
EpsCellData cd(*c);
cells2.swap (cells);
-// Next we have to shift and scale
+ // Next we have to shift and scale
// a bit so that everything still
// arrives in our bounding box of
// 310x310.
};
-// Now we are ready to output...
+ // Now we are ready to output...
for (multiset<DataOut_Old<2>::EpsCellData>::iterator c=cells.begin();
c!=cells.end(); ++c)
{
<< endl;
-///////////////////////////////
+ ///////////////////////////////
// first make up a list of used
// vertices along with their
// coordinates. since dofs have
};
-/////////////////////////////////////
+ /////////////////////////////////////
// now for the cells. this is simpler
// than the above task
if (true)
};
-///////////////////////////////////////
+ ///////////////////////////////////////
// data output.
out << "variable" << endl;
};
-// Cell data is the simplest since the order
+ // Cell data is the simplest since the order
// already is correct
if (true)
{
++n_patches;
-// before we start the loop:
+ // before we start the loop:
// create a quadrature rule that
// actually has the points on this
// patch, and an object that
};
-// now fill in the the data values.
+ // now fill in the the data values.
// note that the required order is
// with highest coordinate running
// fastest, we need to enter each
);
-for (unsigned int i=1; i<values.size(); ++i)
+ for (unsigned int i=1; i<values.size(); ++i)
{
min_value = min (min_value,
*min_element(values[i].begin(),
max_value = min_value+1;
-// now set up the intervals based on
+ // now set up the intervals based on
// the min and max values
intervals.clear ();
// set up one list of intervals
intervals.push_back (intervals[0]);
-// finally fill the intervals
+ // finally fill the intervals
for (unsigned int i=0; i<values.size(); ++i)
for (Vector<number>::const_iterator p=values[i].begin();
p < values[i].end(); ++p)
n_coarsen_fathers);
-// we need counters for
+ // we need counters for
// the 'to_stay_or_refine' cells 'n_sr',
// the 'coarsen_fathers' cells 'n_cf',
// and all the cells where a
};
-const unsigned int dofs_per_cell=dof_handler->get_fe().dofs_per_cell;
+ const unsigned int dofs_per_cell=dof_handler->get_fe().dofs_per_cell;
Vector<number> local_values(dofs_per_cell);
vector<unsigned int> *indexptr;
};
-// append result of this cell
+ // append result of this cell
// to the end of the vector
difference(index) = diff;
};