From: Peter Munch Date: Sun, 8 Dec 2019 19:08:27 +0000 (+0100) Subject: Make mapping available of neighbor during ECL X-Git-Tag: v9.2.0-rc1~749^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=4b53c7c7170161da8494b0b1cc3a04c5064ed155;p=dealii.git Make mapping available of neighbor during ECL --- diff --git a/include/deal.II/matrix_free/fe_evaluation.h b/include/deal.II/matrix_free/fe_evaluation.h index e8bb47f928..c6e2604582 100644 --- a/include/deal.II/matrix_free/fe_evaluation.h +++ b/include/deal.II/matrix_free/fe_evaluation.h @@ -7417,10 +7417,11 @@ FEFaceEvaluationjacobian = &this->matrix_info->get_mapping_info() .face_data_by_cells[this->quad_no] - .jacobians[0][offsets]; - this->normal_x_jacobian = &this->matrix_info->get_mapping_info() - .face_data_by_cells[this->quad_no] - .normals_times_jacobians[0][offsets]; + .jacobians[!this->is_minus_face][offsets]; + this->normal_x_jacobian = + &this->matrix_info->get_mapping_info() + .face_data_by_cells[this->quad_no] + .normals_times_jacobians[!this->is_minus_face][offsets]; # ifdef DEBUG this->dof_values_initialized = false; diff --git a/include/deal.II/matrix_free/mapping_info.templates.h b/include/deal.II/matrix_free/mapping_info.templates.h index 84fdc20559..39be9ccbbe 100644 --- a/include/deal.II/matrix_free/mapping_info.templates.h +++ b/include/deal.II/matrix_free/mapping_info.templates.h @@ -1978,6 +1978,8 @@ namespace internal storage_length * GeometryInfo::faces_per_cell); face_data_by_cells[my_q].jacobians[0].resize_fast( storage_length * GeometryInfo::faces_per_cell); + face_data_by_cells[my_q].jacobians[1].resize_fast( + storage_length * GeometryInfo::faces_per_cell); if (update_flags & update_normal_vectors) face_data_by_cells[my_q].normal_vectors.resize_fast( storage_length * GeometryInfo::faces_per_cell); @@ -1985,6 +1987,10 @@ namespace internal update_flags & update_jacobians) face_data_by_cells[my_q].normals_times_jacobians[0].resize_fast( storage_length * GeometryInfo::faces_per_cell); + if (update_flags & update_normal_vectors && + update_flags & update_jacobians) + face_data_by_cells[my_q].normals_times_jacobians[1].resize_fast( + storage_length * GeometryInfo::faces_per_cell); if (update_flags & update_jacobian_grads) face_data_by_cells[my_q].jacobian_gradients[0].resize_fast( storage_length * GeometryInfo::faces_per_cell); @@ -2002,6 +2008,10 @@ namespace internal fe_face_values(face_data_by_cells.size()); for (unsigned int i = 0; i < fe_face_values.size(); ++i) fe_face_values[i].resize(face_data_by_cells[i].descriptor.size()); + std::vector>>> + fe_face_values_neigh(face_data_by_cells.size()); + for (unsigned int i = 0; i < fe_face_values_neigh.size(); ++i) + fe_face_values_neigh[i].resize(face_data_by_cells[i].descriptor.size()); for (unsigned int cell = 0; cell < cell_type.size(); ++cell) for (unsigned int my_q = 0; my_q < face_data_by_cells.size(); ++my_q) for (unsigned int face = 0; face < GeometryInfo::faces_per_cell; @@ -2014,8 +2024,17 @@ namespace internal dummy_fe, face_data_by_cells[my_q].descriptor[fe_index].quadrature, update_flags)); + if (fe_face_values_neigh[my_q][fe_index].get() == nullptr) + fe_face_values_neigh[my_q][fe_index].reset( + new dealii::FEFaceValues( + mapping, + dummy_fe, + face_data_by_cells[my_q].descriptor[fe_index].quadrature, + update_flags)); dealii::FEFaceValues &fe_val = *fe_face_values[my_q][fe_index]; + dealii::FEFaceValues &fe_val_neigh = + *fe_face_values_neigh[my_q][fe_index]; const unsigned int offset = face_data_by_cells[my_q] .data_index_offsets[cell * GeometryInfo::faces_per_cell + @@ -2029,6 +2048,19 @@ namespace internal cells[cell * vectorization_width + v].second); fe_val.reinit(cell_it, face); + const bool is_local = cell_it->is_locally_owned(); + + if (is_local) + { + auto cell_it_neigh = + cell_it->neighbor_or_periodic_neighbor(face); + fe_val_neigh.reinit(cell_it_neigh, + cell_it->at_boundary() ? + cell_it->periodic_neighbor_face_no( + face) : + cell_it->neighbor_face_no(face)); + } + // copy data for affine data type if (cell_type[cell] <= affine) { @@ -2051,6 +2083,23 @@ namespace internal inv_jac[d][ee]; } } + if (is_local && (update_flags & update_jacobians)) + for (unsigned int q = 0; q < fe_val.n_quadrature_points; + ++q) + { + DerivativeForm<1, dim, dim> inv_jac = + fe_val_neigh.jacobian(q).covariant_form(); + for (unsigned int d = 0; d < dim; ++d) + for (unsigned int e = 0; e < dim; ++e) + { + const unsigned int ee = ExtractFaceHelper:: + reorder_face_derivative_indices(face, + e); + face_data_by_cells[my_q] + .jacobians[1][offset + q][d][e][v] = + inv_jac[d][ee]; + } + } if (update_flags & update_jacobian_grads) { Assert(false, ExcNotImplemented()); @@ -2117,6 +2166,16 @@ namespace internal .normals_times_jacobians[0][offset + q] = face_data_by_cells[my_q].normal_vectors[offset + q] * face_data_by_cells[my_q].jacobians[0][offset + q]; + if (update_flags & update_normal_vectors && + update_flags & update_jacobians) + for (unsigned int q = 0; q < (cell_type[cell] <= affine ? + 1 : + fe_val.n_quadrature_points); + ++q) + face_data_by_cells[my_q] + .normals_times_jacobians[1][offset + q] = + face_data_by_cells[my_q].normal_vectors[offset + q] * + face_data_by_cells[my_q].jacobians[1][offset + q]; } }