From 4928ebe690d844f486a6afd7079654d4c420ff27 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 27 Apr 1998 09:29:25 +0000 Subject: [PATCH] Small changes and doc updates. Changes some casts and use more obvious iterator syntax. git-svn-id: https://svn.dealii.org/trunk@196 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/base/source/parameter_handler.cc | 96 +++++++++---------- deal.II/deal.II/source/dofs/dof_accessor.cc | 1 + .../deal.II/source/dofs/dof_constraints.cc | 96 ++++++++++--------- deal.II/deal.II/source/dofs/dof_handler.cc | 53 +++++----- deal.II/deal.II/source/grid/tria.cc | 41 ++++---- deal.II/deal.II/source/numerics/assembler.cc | 9 +- deal.II/deal.II/source/numerics/base.cc | 59 +++++------- deal.II/deal.II/source/numerics/data_io.cc | 46 ++++++++- 8 files changed, 222 insertions(+), 179 deletions(-) diff --git a/deal.II/base/source/parameter_handler.cc b/deal.II/base/source/parameter_handler.cc index f6f8c01d05..5358032808 100644 --- a/deal.II/base/source/parameter_handler.cc +++ b/deal.II/base/source/parameter_handler.cc @@ -44,7 +44,7 @@ bool ParameterHandler::read_input (istream &input) { // chars is arbitrary but should // be enough char c[10000]; - input.getline ((char *)&c, 10000); + input.getline (&c[0], 10000); line = c; if (!scan_line (line, lineno)) status = false; @@ -114,9 +114,9 @@ void ParameterHandler::clear () { map >::iterator p; for (p=defaults.subsections.begin(); p!=defaults.subsections.end(); ++p) - delete (*p).second; + delete p->second; for (p=changed_entries.subsections.begin(); p!=changed_entries.subsections.end(); ++p) - delete (*p).second; + delete p->second; defaults.subsections.erase (defaults.subsections.begin(), defaults.subsections.end()); changed_entries.subsections.erase (changed_entries.subsections.begin(), @@ -215,11 +215,11 @@ const String & ParameterHandler::get (const String &entry_string) const { map, less >::const_iterator ptr; ptr = pc->entries.find (entry_string); if (ptr != pc->entries.end()) - return (*ptr).second.first; + return ptr->second.first; // not changed ptr = pd->entries.find (entry_string); - return (*ptr).second.first; + return ptr->second.first; }; @@ -307,32 +307,32 @@ void ParameterHandler::print_parameters_section (ostream &out, // first find out the longest entry name unsigned int longest_entry = 0; for (ptr = pd->entries.begin(); ptr != pd->entries.end(); ++ptr) - if ((*ptr).first.length() > longest_entry) - longest_entry = (*ptr).first.length(); + if (ptr->first.length() > longest_entry) + longest_entry = ptr->first.length(); // print entries one by one for (ptr = pd->entries.begin(); ptr != pd->entries.end(); ++ptr) { // check whether this entry is listed // in the Changed tree - if (pc->entries.find((*ptr).first) != pc->entries.end()) + if (pc->entries.find(ptr->first) != pc->entries.end()) switch (style) { case Text: out << setw(indent_level*2) << "" - << (*ptr).first - << setw(longest_entry-(*ptr).first.length()+3) << " = " - << pc->entries[(*ptr).first].first + << ptr->first + << setw(longest_entry-ptr->first.length()+3) << " = " + << pc->entries[ptr->first].first << " <" - << pd->entries[(*ptr).first].first + << pd->entries[ptr->first].first << ">" << endl; break; case LaTeX: - out << "\\item {\\bf " << (*ptr).first << ":} " - << pc->entries[(*ptr).first].first + out << "\\item {\\bf " << ptr->first << ":} " + << pc->entries[ptr->first].first << " ({\\it default:} " - << pd->entries[(*ptr).first].first + << pd->entries[ptr->first].first << ")" << endl; break; @@ -345,13 +345,13 @@ void ParameterHandler::print_parameters_section (ostream &out, { case Text: out << setw(indent_level*2) << "" - << (*ptr).first - << setw(longest_entry-(*ptr).first.length()+2) << "= " - << (*ptr).second.first << endl; + << ptr->first + << setw(longest_entry-ptr->first.length()+2) << "= " + << ptr->second.first << endl; break; case LaTeX: - out << "\\item {\\bf " << (*ptr).first << ":} " - << (*ptr).second.first + out << "\\item {\\bf " << ptr->first << ":} " + << ptr->second.first << endl; break; default: @@ -368,12 +368,12 @@ void ParameterHandler::print_parameters_section (ostream &out, { case Text: out << setw(indent_level*2) << "" - << "subsection " << (*ptrss).first << endl; + << "subsection " << ptrss->first << endl; break; case LaTeX: out << endl << "\\item {\\bf " - << "Subsection " << (*ptrss).first + << "Subsection " << ptrss->first << "}" << endl << "\\begin{itemize}" << endl; @@ -381,7 +381,7 @@ void ParameterHandler::print_parameters_section (ostream &out, default: Assert (false, ExcNotImplemented()); }; - enter_subsection ((*ptrss).first); + enter_subsection (ptrss->first); print_parameters_section (out, style, indent_level+1); leave_subsection (); switch (style) @@ -552,7 +552,7 @@ ParameterHandler::Section* ParameterHandler::get_present_defaults_subsection () const ParameterHandler::Section* ParameterHandler::get_present_defaults_subsection () const { - Section* sec = (Section*)&defaults; // not nice, but needs to be and + Section* sec = const_cast(&defaults); // not nice, but needs to be and // after all: we do not change #sec# vector::const_iterator SecName = subsection_path.begin(); @@ -583,7 +583,7 @@ ParameterHandler::Section* ParameterHandler::get_present_changed_subsection () { const ParameterHandler::Section* ParameterHandler::get_present_changed_subsection () const { - Section* sec = (Section*)&changed_entries; // same as in get_present_default_s... + Section* sec = const_cast(&changed_entries); // same as in get_present_default_s... vector::const_iterator SecName = subsection_path.begin(); while (SecName != subsection_path.end()) @@ -603,7 +603,7 @@ ParameterHandler::Section::~Section () { map >::iterator p; for (p=subsections.begin(); p!=subsections.end(); ++p) - delete (*p).second; + delete p->second; subsections.erase (subsections.begin(), subsections.end()); @@ -690,7 +690,7 @@ void MultipleParameterLoop::init_branches () { if (multiple_choices.size() > 0) for (vector::iterator i=multiple_choices.end()-1; i >= multiple_choices.begin(); --i) - if ((*i).different_values.size() == 1) + if (i->different_values.size() == 1) multiple_choices.erase (i); // finally calculate number of branches @@ -728,18 +728,18 @@ void MultipleParameterLoop::init_branches_section (const ParameterHandler::Secti // whether it is a multiple entry map, less >::const_iterator e; for (e = sec.entries.begin(); e != sec.entries.end(); ++e) - if ((*e).second.first.contains('{')) + if (e->second.first.contains('{')) multiple_choices.push_back (Entry(subsection_path, - (*e).first, - (*e).second.first)); + e->first, + e->second.first)); // transverse subsections map >::const_iterator s; for (s = sec.subsections.begin(); s != sec.subsections.end(); ++s) { - enter_subsection ((*s).first); - init_branches_section (*(*s).second); + enter_subsection (s->first); + init_branches_section (*s->second); leave_subsection (); }; }; @@ -755,53 +755,53 @@ void MultipleParameterLoop::fill_entry_values (const unsigned int run_no) { { // temporarily enter the subsection tree // of this multiple entry - subsection_path.swap ((*choice).subsection_path); + subsection_path.swap (choice->subsection_path); // set entry Section* pd = get_present_defaults_subsection (); - int selection = (run_no/possibilities) % (*choice).different_values.size(); + int selection = (run_no/possibilities) % choice->different_values.size(); String entry_value; - if ((*choice).type == variant) - entry_value = (*choice).different_values[selection]; + if (choice->type == variant) + entry_value = choice->different_values[selection]; else { - if (run_no>=(*choice).different_values.size()) + if (run_no>=choice->different_values.size()) { cerr << "The given array for entry " - << pd->entries[(*choice).entry_name].first + << pd->entries[choice->entry_name].first << " does not conatin enough elements! Taking empty string instead." << endl; entry_value = ""; } else - entry_value = (*choice).different_values[run_no]; + entry_value = choice->different_values[run_no]; }; // check conformance with regex - if (!entry_value.matches (Regex(pd->entries[(*choice).entry_name].second))) + if (!entry_value.matches (Regex(pd->entries[choice->entry_name].second))) { cerr << "In run no. " << run_no+1 << ":" << endl << " The entry value" << endl << " " << entry_value << endl << " for the entry named" << endl - << " " << (*choice).entry_name << endl + << " " << choice->entry_name << endl << " does not match the given regular expression" << endl - << " " << pd->entries[(*choice).entry_name].second << endl + << " " << pd->entries[choice->entry_name].second << endl << " Taking default value" << endl - << " " << pd->entries[(*choice).entry_name].first << endl; + << " " << pd->entries[choice->entry_name].first << endl; // select default instead - entry_value = pd->entries[(*choice).entry_name].first; + entry_value = pd->entries[choice->entry_name].first; }; Section* pc = get_present_changed_subsection (); - pc->entries[(*choice).entry_name] = make_pair(entry_value, String("")); + pc->entries[choice->entry_name] = make_pair(entry_value, String("")); // get out of subsection again - subsection_path.swap ((*choice).subsection_path); + subsection_path.swap (choice->subsection_path); // move ahead if it was a variant entry - if ((*choice).type == variant) - possibilities *= (*choice).different_values.size(); + if (choice->type == variant) + possibilities *= choice->different_values.size(); }; }; diff --git a/deal.II/deal.II/source/dofs/dof_accessor.cc b/deal.II/deal.II/source/dofs/dof_accessor.cc index 19b39cd2c2..1aee0a8d59 100644 --- a/deal.II/deal.II/source/dofs/dof_accessor.cc +++ b/deal.II/deal.II/source/dofs/dof_accessor.cc @@ -5,6 +5,7 @@ #include #include #include +#include #include diff --git a/deal.II/deal.II/source/dofs/dof_constraints.cc b/deal.II/deal.II/source/dofs/dof_constraints.cc index fed53501ff..c418c3cfce 100644 --- a/deal.II/deal.II/source/dofs/dof_constraints.cc +++ b/deal.II/deal.II/source/dofs/dof_constraints.cc @@ -49,12 +49,12 @@ void ConstraintMatrix::add_entry (const unsigned int line, else { for (line_ptr = &lines.back()-1; line_ptr>=lines.begin(); --line_ptr) - if ((*line_ptr).line == line) + if (line_ptr->line == line) break; Assert (false, ExcLineInexistant(line)); }; - (*line_ptr).entries.push_back (make_pair(column,value)); + line_ptr->entries.push_back (make_pair(column,value)); }; @@ -67,7 +67,7 @@ void ConstraintMatrix::close () { vector::iterator line = lines.begin(), endl = lines.end(); for (; line!=endl; ++line) - sort ((*line).entries.begin(), (*line).entries.end()); + sort (line->entries.begin(), line->entries.end()); // sort the lines sort (lines.begin(), lines.end()); @@ -110,7 +110,7 @@ void ConstraintMatrix::condense (const dSMatrixStruct &uncondensed, new_line.push_back (row); else for (unsigned int row=0; row!=n_rows; ++row) - if (row == (*next_constraint).line) + if (row == next_constraint->line) { // this line is constrained new_line.push_back (-1); @@ -152,11 +152,11 @@ void ConstraintMatrix::condense (const dSMatrixStruct &uncondensed, // let c point to the constraint // of this column vector::const_iterator c = lines.begin(); - while ((*c).line != (unsigned int)uncondensed.get_column_numbers()[j]) + while (c->line != (unsigned int)uncondensed.get_column_numbers()[j]) ++c; - for (unsigned int q=0; q!=(*c).entries.size(); ++q) - condensed.add (new_line[row], new_line[(*c).entries[q].first]); + for (unsigned int q=0; q!=c->entries.size(); ++q) + condensed.add (new_line[row], new_line[c->entries[q].first]); } else // line must be distributed @@ -166,8 +166,8 @@ void ConstraintMatrix::condense (const dSMatrixStruct &uncondensed, // for each entry: distribute if (new_line[uncondensed.get_column_numbers()[j]] != -1) // column is not constrained - for (unsigned int q=0; q!=(*next_constraint).entries.size(); ++q) - condensed.add (new_line[(*next_constraint).entries[q].first], + for (unsigned int q=0; q!=next_constraint->entries.size(); ++q) + condensed.add (new_line[next_constraint->entries[q].first], new_line[uncondensed.get_column_numbers()[j]]); else @@ -177,12 +177,12 @@ void ConstraintMatrix::condense (const dSMatrixStruct &uncondensed, // let c point to the constraint // of this column vector::const_iterator c = lines.begin(); - while ((*c).line != (unsigned int)uncondensed.get_column_numbers()[j]) ++c; + while (c->line != (unsigned int)uncondensed.get_column_numbers()[j]) ++c; - for (unsigned int p=0; p!=(*c).entries.size(); ++p) - for (unsigned int q=0; q!=(*next_constraint).entries.size(); ++q) - condensed.add (new_line[(*next_constraint).entries[q].first], - new_line[(*c).entries[p].first]); + for (unsigned int p=0; p!=c->entries.size(); ++p) + for (unsigned int q=0; q!=next_constraint->entries.size(); ++q) + condensed.add (new_line[next_constraint->entries[q].first], + new_line[c->entries[p].first]); }; ++next_constraint; @@ -212,7 +212,7 @@ void ConstraintMatrix::condense (dSMatrixStruct &sparsity) const { distribute.resize (sparsity.n_rows(), -1); for (unsigned int c=0; c(c); int n_rows = sparsity.n_rows(); for (int row=0; rowline) { // this line is constrained new_line.push_back (-1); @@ -345,13 +345,14 @@ void ConstraintMatrix::condense (const dSMatrix &uncondensed, // let c point to the constraint // of this column vector::const_iterator c = lines.begin(); - while ((*c).line != (unsigned int)uncondensed_struct.get_column_numbers()[j]) ++c; + while (c->line != (unsigned int)uncondensed_struct.get_column_numbers()[j]) + ++c; - for (unsigned int q=0; q!=(*c).entries.size(); ++q) + for (unsigned int q=0; q!=c->entries.size(); ++q) // distribute to rows with // appropriate weight - condensed.add (new_line[row], new_line[(*c).entries[q].first], - uncondensed.global_entry(j) * (*c).entries[q].second); + condensed.add (new_line[row], new_line[c->entries[q].first], + uncondensed.global_entry(j) * c->entries[q].second); } else // line must be distributed @@ -361,11 +362,11 @@ void ConstraintMatrix::condense (const dSMatrix &uncondensed, // for each column: distribute if (new_line[uncondensed_struct.get_column_numbers()[j]] != -1) // column is not constrained - for (unsigned int q=0; q!=(*next_constraint).entries.size(); ++q) - condensed.add (new_line[(*next_constraint).entries[q].first], + for (unsigned int q=0; q!=next_constraint->entries.size(); ++q) + condensed.add (new_line[next_constraint->entries[q].first], new_line[uncondensed_struct.get_column_numbers()[j]], uncondensed.global_entry(j) * - (*next_constraint).entries[q].second); + next_constraint->entries[q].second); else // not only this line but @@ -374,15 +375,16 @@ void ConstraintMatrix::condense (const dSMatrix &uncondensed, // let c point to the constraint // of this column vector::const_iterator c = lines.begin(); - while ((*c).line != (unsigned int)uncondensed_struct.get_column_numbers()[j]) ++c; + while (c->line != (unsigned int)uncondensed_struct.get_column_numbers()[j]) + ++c; - for (unsigned int p=0; p!=(*c).entries.size(); ++p) - for (unsigned int q=0; q!=(*next_constraint).entries.size(); ++q) - condensed.add (new_line[(*next_constraint).entries[q].first], - new_line[(*c).entries[p].first], + for (unsigned int p=0; p!=c->entries.size(); ++p) + for (unsigned int q=0; q!=next_constraint->entries.size(); ++q) + condensed.add (new_line[next_constraint->entries[q].first], + new_line[c->entries[p].first], uncondensed.global_entry(j) * - (*next_constraint).entries[q].second * - (*c).entries[p].second); + next_constraint->entries[q].second * + c->entries[p].second); }; ++next_constraint; @@ -412,7 +414,7 @@ void ConstraintMatrix::condense (dSMatrix &uncondensed) const { distribute.resize (sparsity.n_rows(), -1); for (unsigned int c=0; c(c); int n_rows = sparsity.n_rows(); for (int row=0; rowline) { // this line is constrained new_line.push_back (-1); @@ -560,10 +562,10 @@ void ConstraintMatrix::condense (const dVector &uncondensed, else // line must be distributed { - for (unsigned int q=0; q!=(*next_constraint).entries.size(); ++q) - condensed(new_line[(*next_constraint).entries[q].first]) + for (unsigned int q=0; q!=next_constraint->entries.size(); ++q) + condensed(new_line[next_constraint->entries[q].first]) += - uncondensed(row) * (*next_constraint).entries[q].second; + uncondensed(row) * next_constraint->entries[q].second; ++next_constraint; }; @@ -580,13 +582,13 @@ void ConstraintMatrix::condense (dVector &vec) const { vector::const_iterator next_constraint = lines.begin(); for (unsigned int row=0; rowline) // line must be distributed { - for (unsigned int q=0; q!=(*next_constraint).entries.size(); ++q) - vec((*next_constraint).entries[q].first) + for (unsigned int q=0; q!=next_constraint->entries.size(); ++q) + vec(next_constraint->entries[q].first) += - vec(row) * (*next_constraint).entries[q].second; + vec(row) * next_constraint->entries[q].second; // set entry to zero vec(row) = 0.; @@ -624,7 +626,7 @@ void ConstraintMatrix::distribute (const dVector &condensed, old_line.push_back (row); else for (unsigned int row=0; row!=n_rows; ++row) - if (row == (*next_constraint).line) + if (row == next_constraint->line) { // this line is constrained old_line.push_back (-1); @@ -661,9 +663,9 @@ void ConstraintMatrix::distribute (const dVector &condensed, // first set it to zero uncondensed(line) = 0.; // then add the different contributions - for (unsigned int i=0; i<(*next_constraint).entries.size(); ++i) - uncondensed(line) += (condensed(old_line[(*next_constraint).entries[i].first]) * - (*next_constraint).entries[i].second); + for (unsigned int i=0; ientries.size(); ++i) + uncondensed(line) += (condensed(old_line[next_constraint->entries[i].first]) * + next_constraint->entries[i].second); ++next_constraint; }; }; @@ -678,11 +680,11 @@ void ConstraintMatrix::distribute (dVector &vec) const { { // make entry in line next_constraint.line // first set it to zero - vec((*next_constraint).line) = 0.; + vec(next_constraint->line) = 0.; // then add the different contributions - for (unsigned int i=0; i<(*next_constraint).entries.size(); ++i) - vec((*next_constraint).line) += (vec((*next_constraint).entries[i].first) * - (*next_constraint).entries[i].second); + for (unsigned int i=0; ientries.size(); ++i) + vec(next_constraint->line) += (vec(next_constraint->entries[i].first) * + next_constraint->entries[i].second); }; }; diff --git a/deal.II/deal.II/source/dofs/dof_handler.cc b/deal.II/deal.II/source/dofs/dof_handler.cc index d5f2541407..563da357c3 100644 --- a/deal.II/deal.II/source/dofs/dof_handler.cc +++ b/deal.II/deal.II/source/dofs/dof_handler.cc @@ -6,6 +6,7 @@ #include #include #include +#include #include #include #include @@ -631,19 +632,19 @@ const Triangulation & DoFHandler::get_tria () const { template -const FiniteElement & DoFHandler::get_selected_fe () const { +const FiniteElementBase & DoFHandler::get_selected_fe () const { return *selected_fe; }; template -void DoFHandler::distribute_dofs (const FiniteElement &fe) { +void DoFHandler::distribute_dofs (const FiniteElementBase &fe) { Assert (tria->n_levels() > 0, ExcInvalidTriangulation()); - reserve_space (fe); if (selected_fe != 0) delete selected_fe; - selected_fe = new FiniteElement(fe); + selected_fe = new FiniteElementBase(fe); + reserve_space (); // clear user flags because we will // need them @@ -653,7 +654,7 @@ void DoFHandler::distribute_dofs (const FiniteElement &fe) { active_cell_iterator cell = begin_active(), endc = end(); for (; cell != endc; ++cell) - next_free_dof = distribute_dofs_on_cell (cell, fe, next_free_dof); + next_free_dof = distribute_dofs_on_cell (cell, next_free_dof); used_dofs = next_free_dof; }; @@ -661,7 +662,6 @@ void DoFHandler::distribute_dofs (const FiniteElement &fe) { int DoFHandler<1>::distribute_dofs_on_cell (active_cell_iterator &cell, - const FiniteElement<1> &fe, unsigned int next_free_dof) { // distribute dofs of vertices @@ -681,11 +681,13 @@ int DoFHandler<1>::distribute_dofs_on_cell (active_cell_iterator &cell, // copy dofs { if (v==0) - for (unsigned int d=0; dset_vertex_dof_index (0, d, neighbor->vertex_dof_index (1, d)); + for (unsigned int d=0; ddofs_per_vertex; ++d) + cell->set_vertex_dof_index (0, d, + neighbor->vertex_dof_index (1, d)); else - for (unsigned int d=0; dset_vertex_dof_index (1, d, neighbor->vertex_dof_index (0, d)); + for (unsigned int d=0; ddofs_per_vertex; ++d) + cell->set_vertex_dof_index (1, d, + neighbor->vertex_dof_index (0, d)); // next neighbor continue; @@ -693,12 +695,12 @@ int DoFHandler<1>::distribute_dofs_on_cell (active_cell_iterator &cell, }; // otherwise: create dofs newly - for (unsigned int d=0; ddofs_per_vertex; ++d) cell->set_vertex_dof_index (v, d, next_free_dof++); }; // dofs of line - for (unsigned int d=0; ddofs_per_line; ++d) cell->set_dof_index (d, next_free_dof++); // note that this cell has been processed @@ -711,16 +713,15 @@ int DoFHandler<1>::distribute_dofs_on_cell (active_cell_iterator &cell, int DoFHandler<2>::distribute_dofs_on_cell (active_cell_iterator &cell, - const FiniteElement<2> &fe, unsigned int next_free_dof) { - if (fe.dofs_per_vertex > 0) + if (selected_fe->dofs_per_vertex > 0) // number dofs on vertices for (unsigned int vertex=0; vertex<4; ++vertex) // check whether dofs for this // vertex have been distributed // (only check the first dof) if (cell->vertex_dof_index(vertex, 0) == -1) - for (unsigned int d=0; ddofs_per_vertex; ++d) cell->set_vertex_dof_index (vertex, d, next_free_dof++); // for the four sides @@ -729,12 +730,12 @@ int DoFHandler<2>::distribute_dofs_on_cell (active_cell_iterator &cell, line_iterator line = cell->line(side); // distribute dofs if necessary - if ((fe.dofs_per_line > 0) && + if ((selected_fe->dofs_per_line > 0) && // check whether line dof is already // numbered (check only first dof) (line->dof_index(0) == -1)) // if not: distribute dofs - for (unsigned int d=0; ddofs_per_line; ++d) line->set_dof_index (d, next_free_dof++); // note if line is subject to constraints @@ -744,8 +745,8 @@ int DoFHandler<2>::distribute_dofs_on_cell (active_cell_iterator &cell, // dofs of quad - if (fe.dofs_per_quad > 0) - for (unsigned int d=0; ddofs_per_quad > 0) + for (unsigned int d=0; ddofs_per_line; ++d) cell->set_dof_index (d, next_free_dof++); @@ -903,7 +904,7 @@ void DoFHandler::renumber_dofs (const RenumberingMethod method, //// multimap >::iterator i; for (i = dofs_by_coordination.begin(); i!=dofs_by_coordination.end(); ++i) - new_number[(*i).second] = next_free_number++; + new_number[i->second] = next_free_number++; // after that: copy this round's // dofs for the next round @@ -1512,7 +1513,8 @@ void DoFHandler::distribute_cell_to_dof_vector (const dVector &cell_data, -void DoFHandler<1>::reserve_space (const FiniteElement<1> &fe) { +void DoFHandler<1>::reserve_space () { + Assert (selected_fe != 0, ExcNoFESelected()); Assert (tria->n_levels() > 0, ExcInvalidTriangulation()); // delete all levels and set them up // newly, since vectors are @@ -1530,14 +1532,15 @@ void DoFHandler<1>::reserve_space (const FiniteElement<1> &fe) { levels.push_back (new DoFLevel<1>); levels.back()->line_dofs.resize (tria->levels[i]->lines.lines.size() * - fe.dofs_per_line, + selected_fe->dofs_per_line, -1); }; }; -void DoFHandler<2>::reserve_space (const FiniteElement<2> &fe) { +void DoFHandler<2>::reserve_space () { + Assert (selected_fe != 0, ExcNoFESelected()); Assert (tria->n_levels() > 0, ExcInvalidTriangulation()); // delete all levels and set them up @@ -1557,10 +1560,10 @@ void DoFHandler<2>::reserve_space (const FiniteElement<2> &fe) { levels.push_back (new DoFLevel<2>); levels.back()->line_dofs.resize (tria->levels[i]->lines.lines.size() * - fe.dofs_per_line, + selected_fe->dofs_per_line, -1); levels.back()->quad_dofs.resize (tria->levels[i]->quads.quads.size() * - fe.dofs_per_quad, + selected_fe->dofs_per_quad, -1); }; }; diff --git a/deal.II/deal.II/source/grid/tria.cc b/deal.II/deal.II/source/grid/tria.cc index cadf461289..43156902bc 100644 --- a/deal.II/deal.II/source/grid/tria.cc +++ b/deal.II/deal.II/source/grid/tria.cc @@ -116,7 +116,7 @@ void Triangulation<1>::create_triangulation (const vector > &v, if (lines_at_vertex[line->vertex_index(vertex)][0] == line->index()) if (lines_at_vertex[line->vertex_index(vertex)].size() == 2) { - const cell_iterator neighbor ((Triangulation<1>*)this, + const cell_iterator neighbor (const_cast*>(this), 0, // level lines_at_vertex[line->vertex_index(vertex)][1]); line->set_neighbor (vertex, neighbor); @@ -129,7 +129,7 @@ void Triangulation<1>::create_triangulation (const vector > &v, // present line is not first adjacent // one -> first adjacent one is neighbor { - const cell_iterator neighbor ((Triangulation<1>*)this, + const cell_iterator neighbor (const_cast*>(this), 0, // level lines_at_vertex[line->vertex_index(vertex)][0]); line->set_neighbor (vertex, neighbor); @@ -266,8 +266,8 @@ void Triangulation<2>::create_triangulation (const vector > &v, for (i=needed_lines.begin(); i!=needed_lines.end(); i++) { // touch the vertices of this line - ++vertex_touch_count[(*i).first.first]; - ++vertex_touch_count[(*i).first.second]; + ++vertex_touch_count[i->first.first]; + ++vertex_touch_count[i->first.second]; }; // assert minimum touch count is at @@ -290,9 +290,9 @@ void Triangulation<2>::create_triangulation (const vector > &v, map,line_iterator,less > >::iterator i; for (i = needed_lines.begin(); line!=end_line(); ++line, ++i) { - line->set (Line((*i).first.first, (*i).first.second)); + line->set (Line(i->first.first, i->first.second)); line->set_used_flag (); - (*i).second = line; + i->second = line; }; }; @@ -367,8 +367,8 @@ void Triangulation<2>::create_triangulation (const vector > &v, for (; boundary_line!=end_boundary_line; ++boundary_line) { line_iterator line; - pair line_vertices(make_pair((*boundary_line).vertices[0], - (*boundary_line).vertices[1])); + pair line_vertices(make_pair(boundary_line->vertices[0], + boundary_line->vertices[1])); if (needed_lines.find(line_vertices) != needed_lines.end()) // line found in this direction line = needed_lines[line_vertices]; @@ -393,7 +393,7 @@ void Triangulation<2>::create_triangulation (const vector > &v, Assert (line->boundary_indicator() == 0, ExcInteriorLineCantBeBoundary()); - line->set_boundary_indicator ((*boundary_line).material_id); + line->set_boundary_indicator (boundary_line->material_id); }; @@ -580,7 +580,7 @@ void Triangulation::save_refine_flags (ostream &out) const { // 3. magic number 0xabcd out << mn_tria_refine_flags_begin << " " << N << endl; for (unsigned int i=0; i(flags[i]) << " "; out << endl; out << mn_tria_refine_flags_end << endl; @@ -685,7 +685,7 @@ void Triangulation::save_user_flags_line (ostream &out) const { // 3. magic number 0xabcf out << mn_tria_line_user_flags_begin << " " << N << endl; for (unsigned int i=0; i(flags[i]) << " "; out << endl; out << mn_tria_line_user_flags_end << endl; @@ -758,7 +758,7 @@ void Triangulation::save_user_flags_quad (ostream &out) const { // 3. magic number 0xabcf out << mn_tria_quad_user_flags_begin << " " << N << endl; for (unsigned int i=0; i(flags[i]) << " "; out << endl; out << mn_tria_quad_user_flags_end << endl; @@ -1119,7 +1119,7 @@ Triangulation::begin_raw_line (unsigned int level) const { if (levels[level]->lines.lines.size() == 0) return end_line (); - return raw_line_iterator ((Triangulation*)this, + return raw_line_iterator (const_cast*>(this), level, 0); }; @@ -1144,7 +1144,7 @@ Triangulation::begin_raw_quad (unsigned int level) const { if (levels[level]->quads.quads.size() == 0) return end_quad(); - return raw_quad_iterator ((Triangulation*)this, + return raw_quad_iterator (const_cast*>(this), level, 0); }; @@ -1233,7 +1233,7 @@ Triangulation::begin_active_quad (unsigned int level) const { template typename TriaDimensionInfo::raw_line_iterator Triangulation::end_line () const { - return raw_line_iterator ((Triangulation*)this, + return raw_line_iterator (const_cast*>(this), -1, -1); }; @@ -1252,7 +1252,7 @@ Triangulation<1>::end_quad () const { template typename TriaDimensionInfo::raw_quad_iterator Triangulation::end_quad () const { - return raw_quad_iterator ((Triangulation*)this, + return raw_quad_iterator (const_cast*>(this), -1, -1); }; @@ -1266,7 +1266,7 @@ Triangulation::last_raw_line (const unsigned int level) const { Assert (levels[level]->lines.lines.size() != 0, ExcEmptyLevel (level)); - return raw_line_iterator ((Triangulation*)this, + return raw_line_iterator (const_cast*>(this), level, levels[level]->lines.lines.size()-1); }; @@ -1291,7 +1291,7 @@ Triangulation::last_raw_quad (const unsigned int level) const { Assert (levels[level]->quads.quads.size() != 0, ExcEmptyLevel (level)); - return raw_quad_iterator ((Triangulation*)this, + return raw_quad_iterator (const_cast*>(this), level, levels[level]->quads.quads.size()-1); }; @@ -1667,8 +1667,9 @@ unsigned int Triangulation::max_adjacent_cells () const { for (unsigned vertex=0; vertex<(1<vertex_index(vertex)]; - return max ((unsigned int)1<(1<(*(max_element (usage_count.begin(), + usage_count.end())))); }; diff --git a/deal.II/deal.II/source/numerics/assembler.cc b/deal.II/deal.II/source/numerics/assembler.cc index 616ac3b0ab..a88b02e521 100644 --- a/deal.II/deal.II/source/numerics/assembler.cc +++ b/deal.II/deal.II/source/numerics/assembler.cc @@ -3,6 +3,7 @@ #include #include #include +#include #include #include #include @@ -51,7 +52,7 @@ AssemblerData::AssemblerData (const DoFHandler &dof, dVector &rhs_vector, const Quadrature &quadrature, const FiniteElement &fe, - const UpdateFields &update_flags, + const UpdateFlags &update_flags, const Boundary &boundary) : dof(dof), assemble_matrix(assemble_matrix), @@ -98,9 +99,9 @@ Assembler::Assembler (Triangulation *tria, template void Assembler::assemble (const Equation &equation) { // re-init fe values for this cell - fe_values.reinit (Triangulation::cell_iterator (tria, - present_level, - present_index), + fe_values.reinit (DoFHandler::cell_iterator (tria, + present_level, + present_index), fe, boundary); const unsigned int n_dofs = dof_handler->get_selected_fe().total_dofs; diff --git a/deal.II/deal.II/source/numerics/base.cc b/deal.II/deal.II/source/numerics/base.cc index 6e0aff2d48..836b632fd2 100644 --- a/deal.II/deal.II/source/numerics/base.cc +++ b/deal.II/deal.II/source/numerics/base.cc @@ -6,6 +6,7 @@ #include #include #include +#include #include #include "../../../mia/control.h" @@ -76,7 +77,7 @@ template void ProblemBase::assemble (const Equation &equation, const Quadrature &quadrature, const FiniteElement &fe, - const UpdateFields &update_flags, + const UpdateFlags &update_flags, const DirichletBC &dirichlet_bc, const Boundary &boundary) { Assert ((tria!=0) && (dof_handler!=0), ExcNoTriaSelected()); @@ -167,30 +168,21 @@ void ProblemBase::integrate_difference (const Function &exact_sol difference.reinit (tria->n_active_cells()); - UpdateFields update_flags = UpdateFields (update_q_points | - update_jacobians | - update_JxW_values); + UpdateFlags update_flags = UpdateFlags (update_q_points | + update_jacobians | + update_JxW_values); if ((norm==H1_seminorm) || (norm==H1_norm)) - update_flags = UpdateFields (update_flags | update_gradients); + update_flags = UpdateFlags (update_flags | update_gradients); FEValues fe_values(fe, q, update_flags); // loop over all cells - // (we need an iterator on the triangulation for - // fe_values.reinit, but we also need an iterator - // on the dofhandler; we could generate the first - // out of the second or vice versa each time - // we need them, but conversion seems to be slow. - // therefore, we keep two iterators, hoping that - // incrementing is faster than conversion...) DoFHandler::active_cell_iterator cell = dof_handler->begin_active(), endc = dof_handler->end(); - Triangulation::active_cell_iterator tria_cell = dof_handler->get_tria().begin_active(); - - for (unsigned int index=0; cell != endc; ++cell, ++tria_cell, ++index) + for (unsigned int index=0; cell != endc; ++cell, ++index) { double diff=0; // initialize for this cell - fe_values.reinit (tria_cell, fe, boundary); + fe_values.reinit (cell, fe, boundary); switch (norm) { @@ -394,7 +386,7 @@ void ProblemBase::apply_dirichlet_bc (dSMatrix &matrix, make_boundary_value_list (dirichlet_bc, fe, boundary, boundary_values); map::const_iterator dof, endd; - const unsigned int n_dofs = (unsigned int)matrix.m(); + const unsigned int n_dofs = matrix.m(); const dSMatrixStruct &sparsity = matrix.get_sparsity_pattern(); const unsigned int *sparsity_rowstart = sparsity.get_rowstart_indices(); const int *sparsity_colnums = sparsity.get_column_numbers(); @@ -405,9 +397,9 @@ void ProblemBase::apply_dirichlet_bc (dSMatrix &matrix, // set entries of this line // to zero - for (unsigned int j=sparsity_rowstart[(*dof).first]; - jfirst]; + jfirst+1]; ++j) + if (sparsity_colnums[j] != dof->first) // if not main diagonal entry matrix.global_entry(j) = 0.; @@ -425,9 +417,9 @@ void ProblemBase::apply_dirichlet_bc (dSMatrix &matrix, // store the new rhs entry to make // the gauss step more efficient double new_rhs; - if (matrix.diag_element((*dof).first) != 0.0) - new_rhs = right_hand_side((*dof).first) - = (*dof).second * matrix.diag_element((*dof).first); + if (matrix.diag_element(dof->first) != 0.0) + new_rhs = right_hand_side(dof->first) + = dof->second * matrix.diag_element(dof->first); else { double first_diagonal_entry = 1; @@ -438,23 +430,23 @@ void ProblemBase::apply_dirichlet_bc (dSMatrix &matrix, break; }; - matrix.set((*dof).first, (*dof).first, + matrix.set(dof->first, dof->first, first_diagonal_entry); - new_rhs = right_hand_side((*dof).first) - = (*dof).second * first_diagonal_entry; + new_rhs = right_hand_side(dof->first) + = dof->second * first_diagonal_entry; }; // store the only nonzero entry // of this line for the Gauss // elimination step - const double diagonal_entry = matrix.diag_element((*dof).first); + const double diagonal_entry = matrix.diag_element(dof->first); // do the Gauss step for (unsigned int row=0; rowfirst) && + ((signed int)row != dof->first)) // this line has an entry // in the regarding column // but this is not the main @@ -470,7 +462,7 @@ void ProblemBase::apply_dirichlet_bc (dSMatrix &matrix, // preset solution vector - solution((*dof).first) = (*dof).second; + solution(dof->first) = dof->second; }; }; @@ -503,7 +495,6 @@ ProblemBase::make_boundary_value_list (const DirichletBC &dirichlet_ // a Tria-iterator for the fe object DoFHandler::active_face_iterator face = dof_handler->begin_active_face(), endf = dof_handler->end_face(); - Triangulation::active_face_iterator tface = tria->begin_active_face(); DirichletBC::const_iterator function_ptr; @@ -515,7 +506,7 @@ ProblemBase::make_boundary_value_list (const DirichletBC &dirichlet_ vector > dof_locations (face_dofs.size(), Point()); vector dof_values; - for (; face!=endf; ++face, ++tface) + for (; face!=endf; ++face) if ((function_ptr = dirichlet_bc.find(face->boundary_indicator())) != dirichlet_bc.end()) // face is subject to one of the @@ -527,8 +518,8 @@ ProblemBase::make_boundary_value_list (const DirichletBC &dirichlet_ face_dofs.erase (face_dofs.begin(), face_dofs.end()); dof_values.erase (dof_values.begin(), dof_values.end()); face->get_dof_indices (face_dofs); - fe.face_ansatz_points (tface, boundary, dof_locations); - (*function_ptr).second->value_list (dof_locations, dof_values); + fe.get_face_ansatz_points (face, boundary, dof_locations); + function_ptr->second->value_list (dof_locations, dof_values); // enter into list for (unsigned int i=0; i #include #include +#include + #include #include #include @@ -435,6 +437,48 @@ void DataOut::write_gnuplot (ostream &out) const { Assert (dofs->get_selected_fe().dofs_per_vertex==1, ExcIncorrectDofsPerVertex()); + // write preamble + if (true) + { + // block this to have local + // variables destroyed after + // use + time_t time1= time (0); + tm *time = localtime(&time1); + out << "# This file was generated by the deal.II library." << endl + << "# Date = " + << time->tm_year+1900 << "/" + << time->tm_mon+1 << "/" + << time->tm_mday << endl + << "# Time = " + << time->tm_hour << ":" + << setw(2) << time->tm_min << ":" + << setw(2) << time->tm_sec << endl + << "#" << endl + << "# For a description of the UCD format see the AVS Developers guide." + << endl + << "#" << endl + << "# "; + switch (dim) + { + case 1: + out << " "; + break; + case 2: + out << " "; + break; + default: + Assert (false, ExcNotImplemented()); + }; + for (unsigned int i=0; i!=data.size(); ++i) + out << '<' + << data[i].name << ':' << data[i].units + << "> "; + out << endl; + + }; + + DoFHandler::active_cell_iterator cell; DoFHandler::active_cell_iterator endc = dofs->end(); @@ -483,7 +527,7 @@ void DataOut::write_gnuplot (ostream &out) const { break; - case 3: + default: Assert (false, ExcNotImplemented()); }; }; -- 2.39.5