Dear deal.ii users
A question concerning the implementation of a problem involving multiple
degrees of freedom per node. I'm solving a problem in nonlocal crystal
plasticity where the nodal degrees of freedom are the displacement (3
components per node in R^3) and the plastic slip. The number of nodal unknowns
for the plastic slip depends on the type of crystal and is therefore decided
at run-time. For, example an fcc crystal would have 12 slip degrees of freedom
per node.
The resulting system of non-linear equations is most conveniently constructed
using Block structures and I therefore follow the approach adopted in Steps
20-21 and construct the fe space as follows:
fe(FE_Q<deal_II_dimension>(pOrderIn), deal_II_dimension,
FE_Q<deal_II_dimension>(pOrderIn), crystal_in.get_no_slip_systems()),
where pOrderIn is the order of the polynomial and
crystal_in.get_no_slip_systems() returns the number of slip systems.
My questions concerns the use of the FEValuesExtractor class. I set up the
extractors as follows:
const FEValuesExtractors::Vector displacements(0);
const FEValuesExtractors::Vector slips(deal_II_dimension);
This works well but what is the most convenient way to access the values of
the shape functions (and their gradients) for the shape functions related to
the slips?
If I use :
const Tensor < 1, deal_II_dimension> phi_Gamma_i = fe_values[slips].value(i,
q);
I get a 1st order tensor with deal_II_dimension = 3 components, but what does
this imply as I have any number of slip systems? So it does not appear that
this is the way to extract the data.
I've experimented using the "fe.system_to_base_index()" command to determine
if the shape functions correspond to the displacement or the slips. Is this
"safe" when I have renumbered the nodes for the block structure? (i.e. will
the line noted ** actually produce the first non-zero value of the shape
function associated with the slip shape functions ?). My numerical trials
indicate it does not.
for (unsigned int q = 0; q < n_q_points; ++q)
for (unsigned int i = 0; i < dofs_per_cell; ++i) {
// 0 = U, 1 = Gamma interpolation fields
const unsigned int shape_func_i_group =
fe.system_to_base_index(i).first.first;
// the index of the non-zero shape function
const unsigned int shape_func_i_index =
fe.system_to_base_index(i).first.second;
for (unsigned int j = 0; j < dofs_per_cell; ++j) {
// 0 = U, 1 = Gamma interpolation fields
const unsigned int shape_func_j_group =
fe.system_to_base_index(j).first.first;
// the index of the non-zero shape function
const unsigned int shape_func_j_index =
fe.system_to_base_index(j).first.second;
if ( (shape_func_i_group == shape_func_j_group) &&
(shape_func_i_group == 0) ) { /
// U - U block
} else if ( (shape_func_i_group == shape_func_j_group) &&
(shape_func_i_group == 1) ) {
// Gamma - Gamma block
} else if ( (shape_func_i_group == 0) &&
(shape_func_j_group == 1) ) {
// U - Gamma block
double phi_Gamma_j = fe_values.shape_value(j, q); // (**)
} else if ( (shape_func_i_group == 1) &&
(shape_func_j_group == 0) ) {
// Gamma - U block
}
}
}
Any help would be most appreciated.
Regards
Andrew
_______________________________________________
