Daniel, thank you for the explanation. I understand that the components are in space (x/y/z), and (*shape_gradient_ptr++)[d] would be derivative of phi in each direction; is that correct?
(compared to vectors) the implementation of divergence function for the tensors is as following; *if* (snc != -1) { *const* *unsigned* *int* comp = shape_function_data[shape_function].single_nonzero_component_index; *const* dealii::Tensor < 1, *spacedim*> *shape_gradient_ptr = &shape_gradients[snc][0]; *const* TableIndices<2> indices = dealii::Tensor<2,*spacedim* >::unrolled_to_component_indices(comp); *const* *unsigned* *int* ii = indices[0]; *const* *unsigned* *int* jj = indices[1]; *for* (*unsigned* *int* q_point = 0; q_point < n_quadrature_points; ++q_point, ++shape_gradient_ptr) { divergences[q_point][jj] += value * (*shape_gradient_ptr)[ii]; } } *else* { *for* (*unsigned* *int* d = 0; d < *dim***dim*; ++d) *if* (shape_function_data[shape_function].is_nonzero_shape_function_component[d]) { Assert (*false*, ExcNotImplemented()); } } There is only the "snc != 1" case, and it seems that the component now means the tensor components and not the directions in space. Is this function (for tensors) not yet implemented? or is it sth related to non-primitive shape functions (I would appreciate if you could comment on this as well)? Thank you in advance. Best, Metin On Sunday, October 2, 2016 at 1:13:55 PM UTC+2, Daniel Arndt wrote: > > Metin, > > [...] >> For the vectors, the implementation (for the else case) is as following; >> >> *for* (*unsigned* *int* shape_function=0; >> >> shape_function<dofs_per_cell; ++shape_function) >> >> { >> >> *for* (*unsigned* *int* d=0; d<*spacedim*; ++d) >> >> *if* >> (shape_function_data[shape_function].is_nonzero_shape_function_component[d]) >> >> { >> >> *const* dealii::Tensor<1,*spacedim*> >> *shape_gradient_ptr = >> >> &shape_gradients[shape_function_data[shape_function]. >> >> row_index[d]][0]; >> >> *for* (*unsigned* *int* q_point=0; >> q_point<n_quadrature_points; ++q_point) >> >> divergences[q_point] += value * >> (*shape_gradient_ptr++)[d]; >> >> } >> >> } >> >> To me, it's not clear how the "value and shape_gradient_ptr" is used or >> in general how the (e.g.) u_x+v_y+w_z is done.. >> > Let me just comment on this: > The divergence at a quadrature point is decomposed into a sum over all > shape functions and all of their components. > Since the value at a certain quadrature point q_p of a FE function u can > be written as u(q_p)=\sum_i u_i \phi_i(q_p), > it holds \nabla \cdot u(q_p) = \sum_i u_i \nabla \cdot \phi_i(q_p) = > \sum_i u_i \sum_c \partial_c \phi_i^c(q_p) where \phi_i^c is the cth > component of the ith shape function. > This is exactly what is done in the last line. > shape_gradient_ptr is initialized with the gradient of the dth component > of the considered shape_function at the first quadrature point. > By "*shape_gradient_ptr++", the gradient at this quadrature point is > returned and the pointer incremented such that it points to the > gradient at the next quadrature point. > > Best, > Daniel > -- The deal.II project is located at http://www.dealii.org/ For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en --- You received this message because you are subscribed to the Google Groups "deal.II User Group" group. To unsubscribe from this group and stop receiving emails from it, send an email to dealii+unsubscr...@googlegroups.com. For more options, visit https://groups.google.com/d/optout.