Hi Re the issue of projecting tensors stored at quadrature points to nodes:
I'm not sure this is worth considering adding to deal.ii, but i wrote two functions (attached) that project first-order and symmetric second-order tensors from quadrature points to the nodes of the cell. They are the extension that is proposed in step-18 for dealing with problems where history variables are stored at quadrature points. They extend the functionality of FETools::compute_projection_from_quadrature_points_matrix The attached test of the functionality generates data at a quadrature point and then projects to the nodes. As a check (these things are not always obvious) it interpolates the projected nodal values (only the zeroth component actually) back to the quadrature points. It all seems to work fine. Wolfgang, let me know if you think it's worth adding or needs some reworking... Regards Andrew
fe_tools_interpolate_tensors.tgz
Description: Binary data
On 28 Jun 2010, at 4:47 PM, Michael Rapson wrote:
> Hi there Andrew,
>
> I have had a similar problem when I was trying to project the gradient
> of a complex scalar field (implemented as a vector field with two
> components) to the nodes. The way that I solved my problem was to take
> the original projection matrix and copy it twice into a new matrix of
> twice the size. The various elements in the original matrix need to be
> placed into the correct places in the new matrix. This basically
> relies on the fact that the projections for each of the components are
> independent.
>
> My code to distribute the matrix, which seems to be generalized to
> more components, was as follows:
> // find flows from the gradient in the x and y direction
> {
> QGauss<dim> lhs_quadrature(quadrature_order);
> QGauss<dim> rhs_quadrature(quadrature_order);
> FullMatrix<double> transform_matrix;
>
> FETools::compute_projection_from_quadrature_points_matrix
> (FE_Q<dim>(fe_order),
> lhs_quadrature, rhs_quadrature,
> transform_matrix);
> // transform_matrix.print(std::cout);
> unsigned int n_rows = transform_matrix.size(0);
> unsigned int n_columns = transform_matrix.size(1);
> unsigned int n_components = fe.n_components ();
>
>
> FullMatrix<double> double_transform_matrix
> (n_components*n_rows,
> n_components*n_columns);
> for (unsigned int row = 0; row < n_rows; row++)
> {
> for (unsigned int column = 0; column <
> n_columns; column++)
> {
>
> double_transform_matrix(n_components*row, n_components*column) =
> transform_matrix(row, column);
>
> double_transform_matrix(n_components*row + 1, n_components*column
> + 1) = transform_matrix(row, column);
> }
> }
> // double_transform_matrix.print(std::cout);
>
> You'll notice that I am assuming a specific ordering of the components
> and dofs in the FE, i.e. deal.II internals, which is not recommended
> generally. This code worked for 6.1.0 and I presume it will work for
> the more recent versions of the library too, but you should take this
> into account if it is more than just a test.
>
> Secondly I was using the same FE (FE_Q(N)) for all components. I
> presume you'll have fewer pressure components so this code will not
> work straight up.
>
> Cheers hey,
> Michael
>
>
> On Mon, Jun 28, 2010 at 12:00 PM, Andrew McBride
> <[email protected]> wrote:
>> Hi all
>>
>> I'm trying to implement an adaptive strategy for a problem in non-linear
>> solid mechanics where internal variables are stored at the level of the
>> quadrature points. It's in the spirit of the extension of step-18 as
>> proposed in the documentation. The objective is to project data (scalar,
>> vector, tensor) stored at quadrature points between quadrature points on
>> successively refined / coarsened meshes.
>>
>> The first snag I have hit is that I have vector valued shape function for
>> the displacement (in R^3) and the concentration (in R). The linear system is
>> furthermore split into blocks. The underlying shape functions are primitive
>> however.
>>
>> I wish to use the function "
>> FETools::compute_projection_from_quadrature_points_matrix" to perform the
>> first step of projecting the quadrature point data on the old mesh to the
>> nodes of the old mesh. The snag is that "fe" has two components
>> (displacement and concentration resp.) and hence the method fails.
>>
>> i.e. if i use
>> FETools::compute_projection_from_quadrature_points_matrix(fe,
>> quadrature_formula, quadrature_formula, projection_matrix);
>> i get the following error:
>> An error occurred in line <1924> of file
>> </Users/andrewmcbride/lib/real_deal_svn/deal.II/deal.II/source/fe/fe_tools.cc>
>> in function
>> static void
>> dealii::FETools::compute_projection_from_quadrature_points_matrix(const
>> dealii::FiniteElement<dim, spacedim>&, const dealii::Quadrature<dim>&, const
>> dealii::Quadrature<dim>&, dealii::FullMatrix<double>&) [with int dim = 3,
>> int spacedim = 3]
>> The violated condition was:
>> fe.n_components() == 1
>> The name and call sequence of the exception was:
>> ExcNotImplemented()
>> Additional Information:
>> (none)
>>
>> Is there a simple way around this? Can I send a restriction of "fe" into the
>> function " FETools::compute_projection_from_quadrature_points_matrix"?
>>
>> Has anyone else implemented similar functionality into deal.II? It seems to
>> be essential for adaptive problems in (non-linear) solid mechanics. If so,
>> any general suggestions would be appreciated. Perhaps we can formalise this
>> functionality into deal.ii if there is interest.
>>
>> Thanks again
>> Andrew
>> _______________________________________________
>> dealii mailing list http://poisson.dealii.org/mailman/listinfo/dealii
>>
_______________________________________________ dealii mailing list http://poisson.dealii.org/mailman/listinfo/dealii
