Re: [deal.II] Interpolation between FE_Q- and FE_Bernstein elements
On a second test I found out that by resetting the out-vector to 0 I get
the expected output, regardless of direction. Does that make sense?
Am Donnerstag, 30. Mai 2019 22:25:59 UTC+2 schrieb Maxi Miller:
>
> That means that I don't have to create that matrix myself, but rather use
> the matrix I get from FE_Q_*::get_interpolation_matrix()?
> The most interesting bug my code has at the moment is:
> Q(a)->B(b): I get identical vectors, i.e. a == b
> B(b)->Q(a): Vectors are not identical anymore, a != b
> Q(a)->B(b): Even though that worked in the first round, it does not in the
> second, a != b
> I do not understand that behavior yet.
>
> Am Donnerstag, 30. Mai 2019 19:00:02 UTC+2 schrieb Wolfgang Bangerth:
>>
>> On 5/30/19 9:06 AM, 'Maxi Miller' via deal.II User Group wrote:
>> > I wanted to write some functions to convert my solution from using
>> > FE_Q-elements to FE_Bernstein-elements (and back). For that I wrote two
>> > functions, convert_feQ_to_feB() and convert_feB_to_feQ() (as listed in
>> the
>> > attachment). When converting from FE_Q-elements to
>> FE_Bernstein-elements, I
>> > get the same values for input and output, but on the way back I get
>> wrong
>> > values (currently I am just using a vector filled with the value 1),
>> but I can
>> > not find the reason why the conversion back is failing. Did I forget
>> something
>> > in the corresponding function?
>>
>> I have to admit that I don't have the time right now to look at your
>> code, but
>> it's really just an interpolation problem to go from FE_B to FE_Q: you
>> need to
>> evaluate each shape function of the FE_B at the interpolation (=support)
>> points of the FE_Q. This gives you a cell-local matrix of the form
>>B_{ij} = \varphi_{Bernstein,j)(x_{j,Q})
>> (or with indices reversed -- check with the documentation).
>>
>> By the way, the place to implement this is in the
>> FE_Q_Base::get_interpolation_matrix() function for B -> Q (see
>> source/fe/fe_base.cc around line 500), and maybe
>> FE_Q_Bernstein::get_interpolation_matrix() if you come up with a clever
>> way to
>> define the operation. I know that the latter function says that you can't
>> do
>> it, but I suspect that if you at least define it in a way so that the
>> interpolation of a function onto itself is the identity operation, and if
>> the
>> interpolation of a FE_Q function onto FE_Bernstein yields the same,
>> continuous
>> function, then it might still go into that place.
>>
>> Best
>> W.
>>
>> --
>>
>> Wolfgang Bangerth email: bang...@colostate.edu
>> www: http://www.math.colostate.edu/~bangerth/
>>
>>
--
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.
To view this discussion on the web visit
https://groups.google.com/d/msgid/dealii/fe785f12-beee-4f10-827d-5869c34d09ff%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
using namespace dealii;
template
void convert_feB_to_feQ(const DoFHandler &dhB,
const DoFHandler &dhq,
const Vector &in,
Vector &out){
out = 0.;
AssertDimension(in.size(), dhB.n_dofs());
const FiniteElement &feB = dhB.get_fe();
const FiniteElement &feq = dhq.get_fe();
const ComponentMask fe_mask(feB.get_nonzero_components(0).size(), true);
AssertDimension(fe_mask.size(), feB.get_nonzero_components(0).size());
std::vector fe_to_real(fe_mask.size(),
numbers::invalid_unsigned_int);
unsigned int size = 0;
for (unsigned int i = 0; i < fe_mask.size(); ++i)
{
if (fe_mask[i])
fe_to_real[i] = size++;
}
const ComponentMask maskq(dim, true);
FullMatrix transfer(feB.dofs_per_cell, feq.dofs_per_cell);
FullMatrix local_transfer(feq.dofs_per_cell);
const std::vector> &points = feq.get_unit_support_points();
std::vector fe_to_feq(feB.dofs_per_cell,
numbers::invalid_unsigned_int);
unsigned int index = 0;
for (unsigned int i = 0; i < feB.dofs_per_cell; ++i)
if (fe_mask[feB.system_to_component_index(i).first])
fe_to_feq[i] = index++;
Assert(index == feB.dofs_per_cell, ExcNotImplemented());
for (unsigned int j = 0; j < feB.dofs_per_cell; ++j)
{
const unsigned int comp_j = feB.system_to_component_index(j).first;
if (fe_mask[comp_j])
for (unsigned int i = 0; i < points.size(); ++i)
{
if (fe_to_real[comp_j] ==
feB.system_to_component_index(i).first)
local_transfer(i, fe_to_feq[j]
Re: [deal.II] Interpolation between FE_Q- and FE_Bernstein elements
That means that I don't have to create that matrix myself, but rather use
the matrix I get from FE_Q_*::get_interpolation_matrix()?
The most interesting bug my code has at the moment is:
Q(a)->B(b): I get identical vectors, i.e. a == b
B(b)->Q(a): Vectors are not identical anymore, a != b
Q(a)->B(b): Even though that worked in the first round, it does not in the
second, a != b
I do not understand that behavior yet.
Am Donnerstag, 30. Mai 2019 19:00:02 UTC+2 schrieb Wolfgang Bangerth:
>
> On 5/30/19 9:06 AM, 'Maxi Miller' via deal.II User Group wrote:
> > I wanted to write some functions to convert my solution from using
> > FE_Q-elements to FE_Bernstein-elements (and back). For that I wrote two
> > functions, convert_feQ_to_feB() and convert_feB_to_feQ() (as listed in
> the
> > attachment). When converting from FE_Q-elements to
> FE_Bernstein-elements, I
> > get the same values for input and output, but on the way back I get
> wrong
> > values (currently I am just using a vector filled with the value 1), but
> I can
> > not find the reason why the conversion back is failing. Did I forget
> something
> > in the corresponding function?
>
> I have to admit that I don't have the time right now to look at your code,
> but
> it's really just an interpolation problem to go from FE_B to FE_Q: you
> need to
> evaluate each shape function of the FE_B at the interpolation (=support)
> points of the FE_Q. This gives you a cell-local matrix of the form
>B_{ij} = \varphi_{Bernstein,j)(x_{j,Q})
> (or with indices reversed -- check with the documentation).
>
> By the way, the place to implement this is in the
> FE_Q_Base::get_interpolation_matrix() function for B -> Q (see
> source/fe/fe_base.cc around line 500), and maybe
> FE_Q_Bernstein::get_interpolation_matrix() if you come up with a clever
> way to
> define the operation. I know that the latter function says that you can't
> do
> it, but I suspect that if you at least define it in a way so that the
> interpolation of a function onto itself is the identity operation, and if
> the
> interpolation of a FE_Q function onto FE_Bernstein yields the same,
> continuous
> function, then it might still go into that place.
>
> Best
> W.
>
> --
>
> Wolfgang Bangerth email: bang...@colostate.edu
>
> www: http://www.math.colostate.edu/~bangerth/
>
>
--
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.
To view this discussion on the web visit
https://groups.google.com/d/msgid/dealii/b7102277-a66f-4e38-84a8-cd90326dca0a%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.
Re: [deal.II] Interpolation between FE_Q- and FE_Bernstein elements
On 5/30/19 9:06 AM, 'Maxi Miller' via deal.II User Group wrote:
> I wanted to write some functions to convert my solution from using
> FE_Q-elements to FE_Bernstein-elements (and back). For that I wrote two
> functions, convert_feQ_to_feB() and convert_feB_to_feQ() (as listed in the
> attachment). When converting from FE_Q-elements to FE_Bernstein-elements, I
> get the same values for input and output, but on the way back I get wrong
> values (currently I am just using a vector filled with the value 1), but I
> can
> not find the reason why the conversion back is failing. Did I forget
> something
> in the corresponding function?
I have to admit that I don't have the time right now to look at your code, but
it's really just an interpolation problem to go from FE_B to FE_Q: you need to
evaluate each shape function of the FE_B at the interpolation (=support)
points of the FE_Q. This gives you a cell-local matrix of the form
B_{ij} = \varphi_{Bernstein,j)(x_{j,Q})
(or with indices reversed -- check with the documentation).
By the way, the place to implement this is in the
FE_Q_Base::get_interpolation_matrix() function for B -> Q (see
source/fe/fe_base.cc around line 500), and maybe
FE_Q_Bernstein::get_interpolation_matrix() if you come up with a clever way to
define the operation. I know that the latter function says that you can't do
it, but I suspect that if you at least define it in a way so that the
interpolation of a function onto itself is the identity operation, and if the
interpolation of a FE_Q function onto FE_Bernstein yields the same, continuous
function, then it might still go into that place.
Best
W.
--
Wolfgang Bangerth email: bange...@colostate.edu
www: http://www.math.colostate.edu/~bangerth/
--
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.
To view this discussion on the web visit
https://groups.google.com/d/msgid/dealii/c18480e7-aab9-53ab-1c28-53b68f890765%40colostate.edu.
For more options, visit https://groups.google.com/d/optout.
[deal.II] Interpolation between FE_Q- and FE_Bernstein elements
I wanted to write some functions to convert my solution from using
FE_Q-elements to FE_Bernstein-elements (and back). For that I wrote two
functions, convert_feQ_to_feB() and convert_feB_to_feQ() (as listed in the
attachment). When converting from FE_Q-elements to FE_Bernstein-elements, I
get the same values for input and output, but on the way back I get wrong
values (currently I am just using a vector filled with the value 1), but I
can not find the reason why the conversion back is failing. Did I forget
something in the corresponding function?
Thanks!
--
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.
To view this discussion on the web visit
https://groups.google.com/d/msgid/dealii/e288a022-ff5b-4876-8ca3-9807605ef2a3%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
using namespace dealii;
template
void convert_feB_to_feQ(const DoFHandler &dhB,
const DoFHandler &dhq,
const Vector &in,
Vector &out){
AssertDimension(in.size(), dhB.n_dofs());
const FiniteElement &feB = dhB.get_fe();
const FiniteElement &feq = dhq.get_fe();
const ComponentMask fe_mask(feB.get_nonzero_components(0).size(), true);
AssertDimension(fe_mask.size(), feB.get_nonzero_components(0).size());
std::vector fe_to_real(fe_mask.size(),
numbers::invalid_unsigned_int);
unsigned int size = 0;
for (unsigned int i = 0; i < fe_mask.size(); ++i)
{
if (fe_mask[i])
fe_to_real[i] = size++;
}
const ComponentMask maskq(dim, true);
FullMatrix transfer(feq.dofs_per_cell, feB.dofs_per_cell);
FullMatrix local_transfer(feB.dofs_per_cell);
const std::vector> &points = feq.get_unit_support_points();
std::vector fe_to_feq(feB.dofs_per_cell,
numbers::invalid_unsigned_int);
unsigned int index = 0;
for (unsigned int i = 0; i < feB.dofs_per_cell; ++i)
if (fe_mask[feB.system_to_component_index(i).first])
fe_to_feq[i] = index++;
Assert(index == feB.dofs_per_cell, ExcNotImplemented());
for (unsigned int j = 0; j < feB.dofs_per_cell; ++j)
{
const unsigned int comp_j = feB.system_to_component_index(j).first;
if (fe_mask[comp_j])
for (unsigned int i = 0; i < points.size(); ++i)
{
if (fe_to_real[comp_j] ==
feB.system_to_component_index(i).first)
local_transfer(i, fe_to_feq[j]) =
feq.shape_value(j, points[i]);
}
}
local_transfer.invert(local_transfer);
for (unsigned int i = 0; i < feq.dofs_per_cell; ++i)
if (fe_to_feq[i] != numbers::invalid_unsigned_int)
for (unsigned int j = 0; j < feB.dofs_per_cell; ++j)
transfer(i, j) = local_transfer(fe_to_feq[i], j);
VectorTools::interpolate(dhB, dhq, transfer, in, out);
}
template
void convert_feQ_to_feB(const DoFHandler &dhq,
const DoFHandler &dhb,
const Vector &in,
Vector &out){
AssertDimension(in.size(), dhq.n_dofs());
const FiniteElement &feq = dhq.get_fe();
const FiniteElement &feB = dhb.get_fe();
const ComponentMask fe_mask(feq.get_nonzero_components(0).size(), true);
AssertDimension(fe_mask.size(), feq.get_nonzero_components(0).size());
std::vector fe_to_real(fe_mask.size(),
numbers::invalid_unsigned_int);
unsigned int size = 0;
for (unsigned int i = 0; i < fe_mask.size(); ++i)
{
if (fe_mask[i])
fe_to_real[i] = size++;
}
const ComponentMask maskq(dim, true);
FullMatrix transfer(feq.dofs_per_cell, feB.dofs_per_cell);
FullMatrix local_transfer(feB.dofs_per_cell);
const std::vector> &points = feq.get_unit_support_points();
std::vector fe_to_feq(feq.dofs_per_cell,
numbers::invalid_unsigned_int);
unsigned int index = 0;
for (unsigned int i = 0; i < feq.dofs_per_cell; ++i)
if (fe_mask[feq.system_to_component_index(i).first])
fe_to_feq[i] = index++;
Assert(index == feq.dofs_per_cell, ExcNotImplemented());
for (unsigned int j = 0; j < feq.dofs_per_cell; ++j)
{
const unsigned int comp_j = feq.system_to_component_index(j).first;
if (fe_mask[comp_j])
for (unsigned int i = 0; i < points.size(); ++i)
{
if (fe_to_real[comp_j] ==
feq.system_to_component_index(i).first)
local_transfer(i, fe_to_feq[j]) =
feq.shape_value(j, points[i]);
}
}
local_transfer.invert(local_transfer);
for (unsigned int i = 0; i < feq.dofs_per_cell; ++i)
if (fe_to_feq[i] != numbers::invalid_unsigned_int)
for (unsigned int j = 0; j < feB.dofs_per_cell; ++j)
transfer(i, j) = local_transfer(fe_to_feq[i], j);
VectorTools::interpolate(dhq, dhb, transfer, in, out);
}
