From 770dd881af645d694123a67c472f587aec7b906d Mon Sep 17 00:00:00 2001 From: wolf Date: Tue, 17 Sep 2002 19:12:56 +0000 Subject: [PATCH] Join two loops, which is actually necessary, since we need some information from the start of the loop in the second loop. The mail that explains this is (sorry, only german): --------------------------------- der shapes-Test geht bei mir fehl, und der Grund ist dass in einem FESystem die _Werte_ nicht von den Basiselementen kopiert werden. Jetzt habe ich diesen Teil gestern geaendert, so dass das wohl an mir liegt, aber ich weiss nicht, warum es vorher funktioniert hat. Also, wie haben ein Q2xQ3 FESystem. Fuer beide muss update_values nicht in update_each gesetzt werden, weil sich die Funktionswerte ja nicht aendern. Aber: irgendwann muss es halt mal ins FESystem kopiert werden. Wo/wann passiert das? Wenn ich naemlich das erste Mal in FESystem::compute_fill lande (von FESystem::reinit auf der ersten Zelle), dann habe ich hier (das ist der alte Code, der neue ist etwas anders, aber im wesentlichen der selbe): for (unsigned int base_no=0, comp=0; base_no &base_fe=base_element(base_no); typename FiniteElementBase::InternalDataBase & base_fe_data = fe_data.get_fe_data(base_no); const UpdateFlags base_flags(dim_1==dim ? base_fe_data.current_update_flags() : base_fe_data.update_flags); [da FEValues ist hier dim_1==dim] [... und weiter unten beim kopieren:] if (base_flags & update_values) data.shape_values(system_index, point)= base_data.shape_values(k,point); Hier enthaelt base_flags kein update_values mehr, so dass die Werte nicht kopiert werden. Da das die einzige Stelle in fe_system ist, wo data.shape_values gesetzt werden, macht mich das ganze ein bisschen stutzig: da kann doch nur was kopiert werden, wenn die base_flags fuer diese Zelle update_values enthalten. Also, nachgeschaut: base_fe_data.update_once hat update_values, aber base_fe_data.first_cell==false... Ah, jetzt sehe ich auch was passiert: vorher hatten wir nur eine einzige Schleife ueber die Basiselement, und in der haben wir die Flags geholt, das Element aufgefordert sich zu initialisieren, und gleich kopiert. Jetzt haben wir zwei Schleifen: einmal initialisieren, dann kopieren. Ich weiss damit auch, wie ich das Problem loesen kann (allerdings erst nach dem Mittagessen), aber deshalb mal die Frage: ob wir auf dem ersten Element sitzen, ist eigentlich nicht damit abgehakt, dass sich das Basiselement initialisiert hat. Trotzdem machen alle Elemente das am Ende des ersten Aufrufs von fill_fe_*_values. Ist das der richtige Platz? Und: ganz kurios macht das das FESystem, naemlich irgendwo in der Mitte der compute_fill-Funktion (die von den fill_fe_*_values Funktionen aufgerufen wird), und zwar *im Innern eines Blocks, der in if (flags & (update_values | update_gradients)) steht*!! Das will mir eigentlich nicht ganz richtig scheinen... Ich werde nachher mal einen vorlaeufigen Fix einchecken, damit es wenigstens wieder geht, aber wuerde trotzdem gerne Eure Meinung dazu hoeren! ------------------------------------- ...and... ------------------------------------- Zu ersterer: die update_flags werden auf eine absolut undurchschaubare Art und Weise gesetzt, geaendert, von irgendwelchen Funktionen zurueckgegeben, etc. Ich habe bis heute nicht verstanden, welche Funktion ich denn wann aufrufen muesste, und dass Flags wie first_cell waehrend einer laufenden (aeusseren) Funktion beeinflussen, was man von current_flags zurueckbekommt, passt da sehr gut ins Bild. Ich meine, mal ehrlich: ist es logisch, dass wenn ich current_update_flags aufrufe, dann fill_fe_values, dann noch mal current_update_flags, dass das Ergebnis dann ein anderes ist? fill_fe_values ist eigentlich ein typisches Beispiel einer Funktion, die _Daten_, nicht den _Zustand_ eines Objektes aendert/aendern sollte! Ich glaube, die Funktion hat bisher eigentlich nur aus Zufall funktioniert, weil eine lokale Variable oben in einem Block deklariert und definiert wurde, waehrend sie erst unten benutzt wurde. Haette man die Deklaration dahin getan, wo die Variable das erste Mal benutzt wird, dann haette das schon nicht mehr funktioniert..... -------------------------------- git-svn-id: https://svn.dealii.org/trunk@6443 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/source/fe/fe_system.cc | 249 ++++++++++++------------- 1 file changed, 117 insertions(+), 132 deletions(-) diff --git a/deal.II/deal.II/source/fe/fe_system.cc b/deal.II/deal.II/source/fe/fe_system.cc index 6d1dc9f681..76dec94c69 100644 --- a/deal.II/deal.II/source/fe/fe_system.cc +++ b/deal.II/deal.II/source/fe/fe_system.cc @@ -573,6 +573,8 @@ FESystem::compute_fill (const Mapping &mapping, typename Mapping::InternalDataBase &fedata, FEValuesData &data) const { + const unsigned int n_q_points = quadrature.n_quadrature_points; + // convert data object to internal // data for this class. fails with // an exception if that is not @@ -621,9 +623,7 @@ FESystem::compute_fill (const Mapping &mapping, // for the base elements. FEValuesData &base_data=fe_data.get_fe_values_data(base_no); const FiniteElement &base_fe=base_element(base_no); - base_data.initialize(quadrature.n_quadrature_points, - base_fe, - base_update_flags); + base_data.initialize (n_q_points, base_fe, base_update_flags); } } @@ -670,6 +670,33 @@ FESystem::compute_fill (const Mapping &mapping, base_fe_data = fe_data.get_fe_data(base_no); FEValuesData & base_data = fe_data.get_fe_values_data(base_no); + + // Make sure that in the + // case of fill_fe_values + // the data is only copied + // from base_data to data + // if base_data is + // changed. therefore use + // fe_fe_data.current_update_flags() + // + // for the case of + // fill_fe_(sub)face_values + // the data needs to be + // copied from base_data to + // data on each face, + // therefore use + // base_fe_data.update_flags. + // + // Store these flags into + // base_flags before + // calling + // base_fe.fill_fe_([sub]face_)values + // as the latter changes + // the return value of + // base_fe_data.current_update_flags() + const UpdateFlags base_flags(dim_1==dim ? + base_fe_data.current_update_flags() : + base_fe_data.update_flags); if (face_no==invalid_face_number) base_fe.fill_fe_values(mapping, cell, @@ -680,134 +707,92 @@ FESystem::compute_fill (const Mapping &mapping, else base_fe.fill_fe_subface_values(mapping, cell, face_no, sub_no, *face_quadrature, mapping_data, base_fe_data, base_data); - }; - // now data has been generated, - // so copy it. we used to work - // by looping over all base - // elements, then over - // multiplicity, then over the - // shape functions from that - // base element, but that - // requires that we can infer - // the global number of a shape - // function from its number in - // the base element. for that - // we had the - // component_to_system_table. - // - // however, this does of course - // no longer work since we have - // non-primitive elements. so - // we go the other way round: - // loop over all shape - // functions of the composed - // element, and find from where - // to copy the data. to be able - // to cache some data, make - // things a little bit more - // complicated: loop over all - // base elements, then over all - // shape functions of the - // composed element, and only - // treat those shape functions - // that belong to a given base - // element - for (unsigned int base_no=0; base_nodofs_per_cell; - ++system_index) - if (this->system_to_base_table[system_index].first.first == base_no) - { - const FiniteElement & - base_fe = base_element(base_no); - typename FiniteElementBase::InternalDataBase & - base_fe_data = fe_data.get_fe_data(base_no); - FEValuesData & - base_data = fe_data.get_fe_values_data(base_no); - - const unsigned int n_q_points = quadrature.n_quadrature_points; - - // Make sure that in - // the case of - // fill_fe_values the - // data is only copied - // from base_data to - // data if base_data is - // changed. therefore - // use - // fe_fe_data.current_update_flags() - // - // for the case of - // fill_fe_(sub)face_values - // the data needs to be - // copied from - // base_data to data on - // each face, therefore - // use - // base_fe_data.update_flags. - // - // Store these flags into - // base_flags before - // calling - // base_fe.fill_fe_([sub]face_)values - // as the latter changes - // the return value of - // base_fe_data.current_update_flags() - const UpdateFlags base_flags(dim_1==dim ? - base_fe_data.current_update_flags() : - base_fe_data.update_flags); - - const unsigned int k = system_to_base_table[system_index].second; - Assert (kn_nonzero_components(i); - unsigned int in_index = 0; - for (unsigned int i=0; in_nonzero_components(system_index) == - base_fe.n_nonzero_components(k), - ExcInternalError()); - for (unsigned int s=0; sn_nonzero_components(system_index); ++s) - { - if (base_flags & update_values) - for (unsigned int q=0; qdofs_per_cell; + ++system_index) + if (this->system_to_base_table[system_index].first.first == base_no) + { + const unsigned int + base_index = system_to_base_table[system_index].second; + Assert (base_indexn_nonzero_components(i); + unsigned int in_index = 0; + for (unsigned int i=0; in_nonzero_components(system_index) == + base_fe.n_nonzero_components(base_index), + ExcInternalError()); + for (unsigned int s=0; sn_nonzero_components(system_index); ++s) + { + if (base_flags & update_values) + for (unsigned int q=0; q::compute_fill (const Mapping &mapping, unsigned int offset = 0; if (face_no != invalid_face_number) offset = (sub_no == invalid_face_number) - ? face_no * quadrature.n_quadrature_points + ? face_no * n_q_points :(face_no * GeometryInfo::subfaces_per_face - + sub_no) * quadrature.n_quadrature_points; + + sub_no) * n_q_points; compute_2nd (mapping, cell, offset, mapping_data, fe_data, data); } } -- 2.39.5