Dear Daniel,
I wrote this function to get the solution. It is the function that I
showed you in other email about using
get_function_values and VectorTools::point_value.
The inputs mapping, dof_handler, solution are taken from the main
class solving the system.
I mean the solution here is the solution vector I get after solving the system.
I think that there should be a way to scale the entries to get the solution.
Could you please give me the idea to do it?
Best regards,
Pham Ngoc Kien
template<int dim>
void OutputTools<dim>::point_solution(const Mapping<dim> &mapping,
const DoFHandler<dim> &dof_handler,
const Vector<double> &solution,
const Point<dim> &point,
Vector<double> &point_solution) const {
const std::pair<typename DoFHandler<dim>::active_cell_iterator, Point<dim>>
cell_point =
GridTools::find_active_cell_around_point(mapping, dof_handler, point);
const double delta = 1e-12;
Assert(GeometryInfo<dim>::distance_to_unit_cell(cell_point.second) < delta,
ExcInternalError());
const Quadrature<dim>
quadrature(GeometryInfo<dim>::project_to_unit_cell(cell_point.second));
FEValues<dim> fe_values(mapping, dof_handler.get_fe(), quadrature,
update_values | update_gradients |
update_quadrature_points);
fe_values.reinit(cell_point.first);
const FEValuesExtractors::Vector vector_re(0);
const FEValuesExtractors::Vector vector_im(dim);
std::vector<Tensor<1, dim>> temp_solution_re(quadrature.size());
std::vector<Tensor<1, dim>> temp_solution_im(quadrature.size());
//Get e_field at point
fe_values[vector_re].get_function_values(solution, temp_solution_re);
fe_values[vector_im].get_function_values(solution, temp_solution_im);
for (unsigned int j = 0; j < dim; ++j) {
point_solution(j) = temp_solution_re[0][j];
point_solution(j + dim) = temp_solution_im[0][j];
}
//testing the resutls
Vector<double> values (dim*2);
VectorTools::point_value(dof_handler,solution,point,values);
std::cout<<"\n";
std::cout<<"cell: \t"<<cell_point.first<<"\n";
std::cout<<"measuring point:
\t"<<fe_values.get_quadrature_points()[0]<<"\n";
std::cout<<"result with get_function_values real: \t"<<
temp_solution_re[0] <<"\n";
std::cout<<"result with get_function_values imag: \t"<<
temp_solution_im[0] <<"\n";
std::cout<<"result with point_value(): \t"<<values<<"\n";
}
Vào Th 7, 23 thg 3, 2019 vào lúc 05:50 Daniel Arndt <
[email protected]> đã viết:
> Pham,
>
> I have a question that when assembling the right hand side with the
>> formula as:
>> F_i = - \int_Omega \varphi_i(x) Js(x) dx at quadrature point x
>> located on the source edge, and F_i = 0 everywhere else.
>> I actually take the integral over the cell rather than the edge, isn't it?
>> And there are 4 cells in 3D and 2 cells in 2D which contains the edge.
>> Should I assemble the right hand side by the following procedure ?
>> 1. find the cell containing the source edge
>> 2. find the dof in the cell represent for the edge.
>> 3. Loop over the cell quadrature points:
>> if the quadrature point are on the edge then compute:
>> F_i += - Omega \varphi_i(x) Js(x) JxW(x)
>> if it not:
>> F_i += 0
>> The right hand side for other dofs that do not represent for source edge
>> are set to be 0.
>>
>
> Yes, this looks like the standard way but I am not sure what you mean by
> "the dof in the cell represent(ative) for the edge".
> In general, just put Js(x) in and assemble on all the cells. Of course,
> you can skip zero contributions.
>
>
>
>> As the shape functions are scaled by a factor of inverse edge length, not
>> only the cell right hand side but also the cell matrix will change when
>> refining the mesh.
>> I don't know if both of them change, the solution would change with them?
>>
> That should not be the case. Matrix and right-hand side should change
> consistently such that the solution is scaled correctly.
>
>
>> I thought about when refining the mesh, for example, the initial mesh has
>> the source edge with 2 vertices are (0,0,0) and (0,0,1),
>> and the refined mesh has 2 source edges with vertices are (0,0,0),
>> (0,0,0.5) and (0,0,0.5), (0,0,1), respectively.
>> The results get by function
>> VectorTools::point_value(dof_handler,solution,point,values);
>> are not the same for the 2 mesh with my codes.
>>
> How do you obtain "values"? The entries need to be scaled accordingly to
> the scaling of the shape functions for representing the same function.
> If "values" is the solution of some (discretized) linear system, the
> difference in the result of VectorTools::point_value should be small after
> mesh refinement.
>
> 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 [email protected].
> For more options, visit https://groups.google.com/d/optout.
>
--
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 [email protected].
For more options, visit https://groups.google.com/d/optout.
