from:"Phạm Ngọc Kiên"

Re: [deal.II] make hanging node constraint for locally refined mesh

2020-05-26 Thread Phạm Ngọc Kiên

Dear Prof. Wolfgang Bangerth,

Thank you very much.

I will change the codes as your guidance.

Yours sincerely,
Kien

Vào Th 4, 27 thg 5, 2020 vào lúc 12:58 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 5/26/20 7:30 PM, Phạm Ngọc Kiên wrote:
> >
> > Thus, the conforming FE_Nedelec and FE_NedelecSZ can only work with the
> > conforming mesh. The non-conforming mesh, which has hanging nodes, does
> not
> > work with these elements.
>
> The FE_NedelecSZ element doesnt, but FE_Nedelec does work with hanging
> nodes.
>
> 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/44ffc21e-c3e2-4f48-a4ec-1537f37e8b42%40colostate.edu
> .
>

-- 
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/CAAo%2BSZcyBae4j_-0ahdsQSrVCbTCrrVOHpixTGiraB8_zONUjg%40mail.gmail.com.

Re: [deal.II] make hanging node constraint for locally refined mesh

2020-05-26 Thread Phạm Ngọc Kiên

Dear Prof. Wolfgang Bangerth,

I used FE_NedelecSZ element.

Thus, the conforming FE_Nedelec and FE_NedelecSZ can only work with the
conforming mesh. The non-conforming mesh, which has hanging nodes, does not
work with these elements.

Thank you very much.

Yours sincerely,
Kien

Vào Th 4, 27 thg 5, 2020 vào lúc 05:42 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 5/26/20 1:34 AM, Phạm Ngọc Kiên wrote:
> > /*  const dealii::FullMatrix& dealii::FiniteElement<,
> >  >::constraints(const dealii::internal::SubfaceCase&)
> const
> > [with int dim = 3; int spacedim = 3]
> > The violated condition was:
> >  (this->dofs_per_face == 0) || (interface_constraints.m() != 0)
> > Additional information:
> >  The finite element for which you try to obtain hanging node
> constraints
> > does not appear to implement them.*/
> > /*
> > */
> > Does this error comes from the AffineConstraints when I try to implement
> the
> > DoFTools::make_hanging_node_constraints?
> >
> > Could you please help me to overcome this error?
>
> Well, the error message says pretty clearly what is happening: The element
> you
> use does not provide hanging node information. What element are you using?
>
> 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/fc63875b-df35-eb80-02cd-c503ca81d82e%40colostate.edu
> .
>

-- 
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/CAAo%2BSZfXpVTE1NW1gMpbsSsng08KXFaJfHqO%2B%2BUVDiMoqFjGXA%40mail.gmail.com.

[deal.II] make hanging node constraint for locally refined mesh

2020-05-26 Thread Phạm Ngọc Kiên

Dear all,
When I try to create a locally refined mesh by the following codes:











* GridGenerator::hyper_cube(triangulation, 0, 1);Point 
center (0.5, 0.5, 0.5);triangulation.refine_global(5);
for (auto cell : triangulation.active_cell_iterators()) {
if (cell->is_locally_owned()) {if 
(cell->point_inside(center)){
cell->set_refine_flag();}}
}triangulation.execute_coarsening_and_refinement();*

The library return this error 




*  const dealii::FullMatrix& dealii::FiniteElement<, 
 >::constraints(const dealii::internal::SubfaceCase&) const 
[with int dim = 3; int spacedim = 3]The violated condition was: 
(this->dofs_per_face == 0) || (interface_constraints.m() != 0)Additional 
information: The finite element for which you try to obtain hanging 
node constraints does not appear to implement them.*

Does this error comes from the AffineConstraints when I try to implement 
the DoFTools::make_hanging_node_constraints?

Could you please help me to overcome this error?

Thank you very much.

Best regards,
Kien

-- 
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/fdcc4fd8-761d-4cc4-a954-8631e0cf117d%40googlegroups.com.

Re: [deal.II] Loading a .msh mesh of big size into parallel::distributed::triangulation using p4est

2019-10-19 Thread Phạm Ngọc Kiên

Hi praveen,
If you run with 1 process(rank = 0) with ./step-40, it's fine.
I want to run with MPI on multiple processors with mpirun -np 4 ./step-40.
In this case the trouble comes.
Thank you very much.
Best
Kien

Vào Th 7, 19 thg 10, 2019 vào lúc 13:45 Praveen C  đã
viết:

> Just an observation running your code with deal.II-9.1.1
>
> In debug mode, it runs at 65 sec.
>
> In release mode, I get segmentation error.
>
> Best
> praveen
>
> $ make
> *Scanning dependencies of target step-40*
> [ 50%] Building CXX object CMakeFiles/step-40.dir/step-40.cc.o
> [100%] *Linking CXX executable step-40*
> [100%] Built target step-40
> $ ./step-40
>
>
> +-+++
> | Total wallclock time elapsed since start|  64.6s ||
> | |||
> | Section | no. calls |  wall time | % of total |
> +-+---+++
> | load grid   | 1 |  64.6s | 1e+02% |
> +-+---+++
>
>
>
>
> +-+++
> | Total wallclock time elapsed since start|  6.34e-06s ||
> | |||
> | Section | no. calls |  wall time | % of total |
> +-+---+++
> +-+---+++
>
> $ make release
> *Scanning dependencies of target release*
> [100%] *Switching CMAKE_BUILD_TYPE to Release*
> -- Autopilot invoked
> -- Run   $ make info  to print a detailed help message
> -- Configuring done
> -- Generating done
> -- Build files have been written to: /Users/praveen/Downloads/loading-test
> ***
> *** Switched to Release mode. Now recompile with:  $ make
> ***
> [100%] Built target release
> $ make
> *Scanning dependencies of target step-40*
> [ 50%] Building CXX object CMakeFiles/step-40.dir/step-40.cc.o
> [100%] *Linking CXX executable step-40*
> [100%] Built target step-40
> $ ./step-40
> Segmentation fault: 11
>
> On 19-Oct-2019, at 4:54 AM, Phạm Ngọc Kiên 
> wrote:
>
> Hi colleagues,
> I have a problem with loading mesh of size more than 100,000 cells using
> MPI
> I am writing to send you my test for loading grid with msh file into a
> parallel::distributed::triangulation.
> The time for loading the grid is extremely huge. It takes most of the
> running time of my code.
> I have tried both on my PC and on our lab's cluster, and it took about 400
> seconds for the msh file containing 131936 cells.
> In case of more than 200,000 cells, I cannot load the mesh.
>
> I also tested that if we load the mesh to Triangulation
> triangulation, the time is much less than that of using
> parallel::distributed::triangulation.
> However, with Triangulation triangulation, I don't think that it is
> suitable for run a large scale problem.
>
>
> --
> 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/C9C81EF1-0A8D-4D7D-8297-88E46465556E%40gmail.com
> <https://groups.google.com/d/msgid/dealii/C9C81EF1-0A8D-4D7D-8297-88E46465556E%40gmail.com?utm_medium=email_source=footer>
> .
>

-- 
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/CAAo%2BSZe%2BaK5O59J824Qzqm_xndSnD2qKFCnYV9unPB8J4StMwg%40mail.gmail.com.

Re: [deal.II] How to get the coordinates of dof with FE_Raviart_Thomas?

2019-09-23 Thread Phạm Ngọc Kiên

Dear Prof. Wolfgang Bangerth,
The function fe.has_support_points() returns false. This means that the
FE_Raviart_Thomas does not have support points.
I think generalized_support_points are used for this type of finite element.

For FE_RaviartThomas
(0),
are the dofs defined at the centers of the faces of a cell?
Is the order of the dofs the same as that of the faces. In that case, I
think I can loop over the faces to evaluate the shape function at a point
on each surface for the purpose of finding the direction.

Thank you very much for your guidance.

Best regards,
Kien


Vào Th 3, 24 thg 9, 2019 vào lúc 12:02 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

>
> Pham,
>
> > With FE_Raviart_Thomas in 3D, the dofs are defined normal to the faces
> of each
> > cell.
> > For example, if we use FE_RaviartThomas
> > <
> https://www.dealii.org/current/doxygen/deal.II/classFE__RaviartThomas.html#adc8e64d0cdb38af6e29dc9ac0279f121>(0),
>
> > then in a cell we have 6 dofs defined on the 6 faces, and there is no
> dofs on
> > the lines.
> > Do we have any way to get the dof coordinates and the direction of the
> shape
> > function of the dof?
>
> You can get the coordinates using DoFTools::map_dofs_to_support_points().
> As
> for the direction: That's something that the RT element does internally
> using
> a convention that takes into account the two adjacent cells. I don't think
> there is a way to determine this direction outside the element short of
> evaluating the shape function at a point on the face.
>
> 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/4457bd19-2f78-fb99-ada7-778b6e921e7e%40colostate.edu
> .
>

-- 
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/CAAo%2BSZfT98hCAOXdpUcAJ%2Bxjr1Bt17AXzY8rJiEj4cb8s7XTAQ%40mail.gmail.com.

[deal.II] How to get the coordinates of dof with FE_Raviart_Thomas?

2019-09-23 Thread Phạm Ngọc Kiên

Hi colleagues,
With FE_Raviart_Thomas in 3D, the dofs are defined normal to the faces of 
each cell.
For example, if we use  FE_RaviartThomas 
(0),
 
then in a cell we have 6 dofs defined on the 6 faces, and there is no dofs 
on the lines.
Do we have any way to get the dof coordinates and the direction of the 
shape function of the dof?
I would like to thank you very much in advance.
Best regards,
Kien

-- 
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/80fb08ab-e39b-4c09-878a-d4cdb5e26562%40googlegroups.com.

Re: [deal.II] How to set material id with MPI

2019-08-29 Thread Phạm Ngọc Kiên

Dear all,
I think we have two ways to do this.
The first one is the way Prof. Wolfgang Bangerth suggested.
The second one is to load the grid in a Triangulation in all processor,
then we set the material id before copying
parallel::distributed::Triangulation from the Triangulation.

When I run the codes in my computer, it takes  a lot of time for p4est to
load the grid.
The loading grid  step is more time consuming than solving the system of
equations with a mesh containing about 100,000 cells.

I would like to thank you very much for your help.
Best,
Kien

Vào Th 7, 24 thg 8, 2019 vào lúc 01:29 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 8/22/19 11:58 PM, Phạm Ngọc Kiên wrote:
> >
> > I have a question for parallel::distributed::Triangulation
> > When 2 cells share 1 edge, but they are living in 2 different MPI
> processes,
> > how can I choose only 1 cell containing the common edge from them.
>
> Is your goal to make sure that only one of the two processors does some
> work
> on these edges? If that's the case, then you need a "tie breaker" -- for
> example, if the subdomain id of a locally owned cell is lower than the
> subdomain of a neighboring ghost cell, then the current processor does the
> work. If the locally owned cell's subdomain id is larger, then the
> neighboring
> processor is in charge of the edge.
>
> 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/a8ecb0ed-359d-fb05-c2b2-e3da8c9218be%40colostate.edu
> .
>

-- 
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/CAAo%2BSZfwoztYKLyFvwRJnwCUKEUf3jM7oF2Fd%2B35F2E4nBRWJw%40mail.gmail.com.

Re: [deal.II] How to set material id with MPI

2019-08-22 Thread Phạm Ngọc Kiên

Hi colleagues,
I have a question for parallel::distributed::Triangulation
When 2 cells share 1 edge, but they are living in 2 different MPI
processes, how can I choose only 1 cell containing the common edge from
them.
I think I have to set material id for the cell in the first process, and
then tell the other one do not set it.
However, I don't know how to send the cell iterator in MPI communication.

Could you please help me to address this issue?
Thank you very much.

Best regards,
Kien



Vào Th 5, 18 thg 7, 2019 vào lúc 17:47 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 7/17/19 7:46 PM, Phạm Ngọc Kiên wrote:
> > I am trying to write codes to find a subset of cells that I want to set
> their
> > material id.
> > The codes run well with 1 processor.
> > However, when testing with more than 1 processor, the codes did wrong
> things.
> > This is because each processor only owns a subset of cells with
> distributed
> > triangulation.
> > Do we have a way to address this issue in deal.II?
>
> In addition to Daniel's questions, take a look at the documentation of the
> parallel::distributed::Triangulation class documentation. It talks about
> similar issues with boundary ids. I would imagine that setting material
> ids
> poses similar challenges, and has similar solutions.
>
> 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/9afe5b13-a5db-ba22-81d1-acd9eb952d7d%40colostate.edu
> .
> 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 dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/CAAo%2BSZfi2S1S2OodmkS8Y-mbnVoSxq361WY9m_r5hFsxdysA6w%40mail.gmail.com.

Re: [deal.II] How to sum up several LinearAlgebraPETSc::MPI::Vector completely_distributed_solution ?

2019-08-02 Thread Phạm Ngọc Kiên

Dear Prof. Wolfgang Bangerth,
I have tried with the sum of the right hand side, but as I think, it did
not work for my problem.
I will try with your guidance.
Thank you very much.
Best regards,
Kien

Vào Th 7, 3 thg 8, 2019 vào lúc 09:42 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 8/1/19 9:24 PM, Phạm Ngọc Kiên wrote:
> >
> > However, I do not know if I can sum up all the
> completely_distributed_solution
> > from the different right-hand-sides in order to get a vector of solution.
> > Could you please tell me how to do this?
>
> You could just sum them up:
>
>AppropriateType sum_of_solutions (...);
>for (unsigned int i=0; i<...; ++i)
>  sum_of_solutions += solution[i];
>
> Of course, the sum of the solutions equals the solution of a single linear
> system with the sum of the right hand sides and the same matrix -- which
> would
> by substantially cheaper to compute.
>
> 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/2871cd63-7236-a03d-cb02-8dc49a98a8db%40colostate.edu
> .
>

-- 
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/CAAo%2BSZeZDo0dLUfG-xtvi1XdRi-BkJb5N0wRCowK52ZrUAa9EA%40mail.gmail.com.

[deal.II] How to sum up several LinearAlgebraPETSc::MPI::Vector completely_distributed_solution ?

2019-08-01 Thread Phạm Ngọc Kiên

Dear colleagues,

I am trying to solve my problem with different right_hand_sides.

For each right hand side, we can solve by the codes:
LinearAlgebraPETSc::MPI::Vector 
completely_distributed_solution(locally_owned_dofs, mpi_communicator);
SolverControl solver_control;
PETScWrappers::SparseDirectMUMPS solver(solver_control, 
mpi_communicator);
solver.solve(system_matrix, completely_distributed_solution, 
system_rhs);
constraints.distribute(completely_distributed_solution);
locally_relevant_solution = completely_distributed_solution;

I think I can write the codes to assemble the different right hand sides to 
solve the system with the same system matrix.
This give me the different solutions.

However, I do not know if I can sum up all the completely_distributed_solution 
from the different right-hand-sides in order to get a vector of solution.
Could you please tell me how to do this? 

I would like to thank you very much for your support.

Best regards,
Kien

-- 
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/d2272747-cbc2-475a-85bc-7fb02f9cce57%40googlegroups.com.

[deal.II] How to set material id with MPI

2019-07-17 Thread Phạm Ngọc Kiên

Hi colleagues,
I am trying to write codes to find a subset of cells that I want to set 
their material id.
The codes run well with 1 processor.
However, when testing with more than 1 processor, the codes did wrong 
things.
This is because each processor only owns a subset of cells with distributed 
triangulation.
Do we have a way to address this issue in deal.II?
Thank you very much.

Best regards,
Kien

-- 
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/2e200ba2-cc84-48ff-b776-c6b1bf144b2e%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] Re: How can we get the factorized system matrix when using SparseDirectMUMPS?

2019-06-28 Thread Phạm Ngọc Kiên

Dear Sebastian,
Thank you very much.
I will try the way you mentioned.
Best regards,
Kien

Vào Th 6, 28 thg 6, 2019 vào lúc 16:31 'Sebastian Stark' via deal.II User
Group  đã viết:

> Kien,
>
>
>> I looked through the documentation in Deal.II before.
>> I think it could be done, too. It would be easier if there is an example
>> to use the PETSc parameters for MATSOLVERMUMPS package in deal.II.
>> Currently, I don't know where to put it in my codes.
>>
>
> As Bruno has said this should work using the command line. When you run
> the program you could e.g. append some option like -mat_mumps_icntl_4 3 to
> set ICNTL(4) to 3. Not sure, but maybe you could do some trickery to modify
> argc and argv if you want to hard code (though this is certainly not a nice
> solution). Anyway, the SparseDirectMUMPS class is rather compact - so if
> you really want to customize, it shouldn't be too hard to write your own
> class using SparseDirectMUMPS as a sample. But it really depends on what
> you want to achieve. If you only want to select the factorization method,
> you can just use set_symmetric_mode(false/true) to switch between "LU" and
> "Cholesky" (which actually is LDL^T).
>
> Regards, Sebastian
>
> --
> 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/13a929d2-45d0-4388-b441-55e149ac9e6e%40googlegroups.com
> 
> .
> 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 dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/CAAo%2BSZcV1F%2B309-LjjeZJ%3DEh%3Dikafx7zdg%3D1JbkTV7VWL6tyrQ%40mail.gmail.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] Re: How can we get the factorized system matrix when using SparseDirectMUMPS?

2019-06-27 Thread Phạm Ngọc Kiên

Bruno,
Thank you very much.
I would check the time for the solve() function with different right hand
side.
If the factorization only be done the first time it is called, the time for
solving the next right hand side should be lesser.

I looked through the documentation in Deal.II before.
I think it could be done, too. It would be easier if there is an example to
use the PETSc parameters for MATSOLVERMUMPS package in deal.II.
Currently, I don't know where to put it in my codes.

Best regards,
Kien

Vào Th 5, 27 thg 6, 2019 vào lúc 21:42 Bruno Turcksin <
bruno.turck...@gmail.com> đã viết:

> Kien,
>
> On Wednesday, June 26, 2019 at 10:20:38 PM UTC-4, Phạm Ngọc Kiên wrote:
>>
>> Hi colleagues,
>> I am trying to build codes using  PETScWrappers::SparseDirectMUMPS
>> solver.
>> Step-62 shows a very good example for using it.
>> However, I really want to know if  we can choose the factorization method
>> in MUMPS by deal.II codes.
>>
> According to the documentation (see here
> <https://www.dealii.org/developer/doxygen/deal.II/classPETScWrappers_1_1SparseDirectMUMPS.html>),
> you can. I don't use PETSc but I think this done using the command line
>
>
>>
>> The second question is: do we have a way to get the factorized system
>> matrix when using this direct solver?
>> Because I want to reuse the factorized system matrix to solve the problem
>> again with other right hand side.
>>
> Just call solve() with a different right-hand-side. The factorization
> should only be done the first time the function is called.
>
> Best,
>
> Bruno
>
> --
> 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/24c660f9-1450-4c1f-b0f7-d9c0a629fde8%40googlegroups.com
> <https://groups.google.com/d/msgid/dealii/24c660f9-1450-4c1f-b0f7-d9c0a629fde8%40googlegroups.com?utm_medium=email_source=footer>
> .
> 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 dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/CAAo%2BSZep%2BHFiY_5St3U-MV%3DRVYK741xJ0QwJWOzG5xy3d-iULQ%40mail.gmail.com.
For more options, visit https://groups.google.com/d/optout.

[deal.II] How can we get the factorized system matrix when using SparseDirectMUMPS?

2019-06-26 Thread Phạm Ngọc Kiên

Hi colleagues,
I am trying to build codes using  PETScWrappers::SparseDirectMUMPS solver.
Step-62 shows a very good example for using it.
However, I really want to know if  we can choose the factorization method 
in MUMPS by deal.II codes.

The second question is: do we have a way to get the factorized system 
matrix when using this direct solver?
Because I want to reuse the factorized system matrix to solve the problem 
again with other right hand side.

I would like to thank you very much for your support.

Best regards,
Kien

-- 
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/1cc4a1aa-078f-4c19-a8f6-5998cbff9819%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] The shape function on physical cell in hexahedral is not parallel with the edge (FeNedelecSZ and FENedelec)

2019-05-07 Thread Phạm Ngọc Kiên

Hi all,
It was a long time since the last time I posted my test for shape function
in 3D here.
This time I send you the results when testing with a more simple model in
2D, and the shape functions are still wrong with non-rectangle mesh.
The codes print results in this form:
dof i_th: 7
edge number: 3
edge vertices:
1 -1 and 0.5 -1
vector from first to second vertex of the edge: -0.5 0
shape functions on real cell:
0 0
0 0
0 0
1.2 1.29551
1.25 0.81524
1.3 0.0147947
2 2.15918
2 1.30438
2 0.022192

Here I use QGaussLobatto to get a certain shape function on the edge of a
cell.
It takes me so many times for testing where my codes return wrong solutions.
I think, in deal.II 9.0, the shape functions, defined in the two classes
FE_Nedelec and FE_NedelecSZ, are only work with structured grid (i.e. cube
in 3D and rectangle in 2D).
I hope it worth mentioning to fix this problem as we all want to carry out
finite element method with unstructured grid.

Best,
Pham Ngoc Kien

Vào Th 3, 16 thg 4, 2019 vào lúc 14:43 Phạm Ngọc Kiên <
ngockien.lea...@gmail.com> đã viết:

> Hi,
> I have written this code for testing purpose.
> In case of a rotated cube, the shape functions seemed to be good as they
> were parallel to the edge, although there was some small round-off errors
> at the value 0.
> However, when I tested with the grid which was loaded from GMSH, the shape
> functions failed.
> In the attachment, I showed the information for the first active cell
> only. It is easy to get information of other cells with the codes.
>
> I think that there exists something wrong with the cell that are not a
> cube, but a hexahedron when initializing the shape functions (fe_values).
>
> I think it would be easier for you to help me with the attachment below.
> Thank you very much.
> Best regards.
> Pham Ngoc Kien
>
> Vào Th 2, 15 thg 4, 2019 vào lúc 21:48 Wolfgang Bangerth <
> bange...@colostate.edu> đã viết:
>
>> On 4/11/19 9:17 PM, Phạm Ngọc Kiên wrote:
>> > Testing for an edge whose global vertices located from (0,0,0) to
>> (0,0,1) in
>> > real coordinates.
>> > With a cube I get the shape function vectors at the dof related to the
>> edge,
>> > for examples, (0,0,0), (0,0,-0.25), (0,0,-0.5), (0,0,-1), which are
>> parallel
>> > to the edge.
>> > However, if I create a mesh with GMSH and check the shape function
>> vectors,
>> > some of them are: (-0.0868383, -0.428, 1), (-0.129933, -0.402772, 1),
>> > (-0.052103 -0.2568 0.6),..,
>> > which are not parallel to the edge.
>> >
>> > In my limited opinion, I will get the wrong solution as the shape
>> functions
>> > created by my code are not parallel to the edge in my model.
>> > Could you please tell me why this happen and how to fix it?
>>
>> Good question. This does sound wrong. To figure out what is happening
>> exactly,
>> we typically construct small test cases that try to illustrate one
>> particular
>> issue. In your case, I would try to create a mesh with exactly one cell
>> --
>> either the one you have above, or maybe even simpler a cube that is just
>> rotated. Then create a DoFHandler on it and output the values of the
>> shape
>> functions. This way, you really have only one thing that could go wrong,
>> and
>> it's easy to demonstrate what you see and how it differs from
>> expectations.
>>
>> Do you think you could come up with a small program that does this?
>>
>> Best
>>   WB
>>
>> --
>> 
>> 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.
>> 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 dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.g

Re: [deal.II] Re: Process finished with exit code 9 ?

2019-05-03 Thread Phạm Ngọc Kiên

Thank you very much.
I will take by myself. I don't think that I am running out of memory
because my computer can solve larger problem previously.
Or somewhere I get the memory leakage.
Best
Pham Ngoc Kien

Vào Th 5, 2 thg 5, 2019 vào lúc 19:40  đã viết:

> Hi,
>
> Is it UMFPACK that returns the exit code? If not it could be that you are
> running out of memory see here
> <https://stackoverflow.com/questions/40888164/c-program-crashes-with-exit-code-9-sigkill>
> .
>
> Best,
>
> Bruno
>
> On Thursday, May 2, 2019 at 2:59:31 AM UTC-4, Phạm Ngọc Kiên wrote:
>>
>> Hi all,
>> When I am trying to solve my problem using the direct solver UMFPACK, my
>> computer can not finish it and return a line:
>> *Process finished with exit code 9*
>> However, the codes did run on my computer before.
>> Could you please tell me the reason?
>> Thank you very much.
>> Pham Ngoc Kien
>>
> --
> 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.
>

-- 
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.

[deal.II] Process finished with exit code 9 ?

2019-05-02 Thread Phạm Ngọc Kiên

Hi all,
When I am trying to solve my problem using the direct solver UMFPACK, my 
computer can not finish it and return a line:
*Process finished with exit code 9*
However, the codes did run on my computer before.
Could you please tell me the reason?
Thank you very much.
Pham Ngoc Kien

-- 
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.

Re: [deal.II] The shape function on physical cell in hexahedral is not parallel with the edge (FeNedelecSZ and FENedelec)

2019-04-15 Thread Phạm Ngọc Kiên

Hi,
I have written this code for testing purpose.
In case of a rotated cube, the shape functions seemed to be good as they
were parallel to the edge, although there was some small round-off errors
at the value 0.
However, when I tested with the grid which was loaded from GMSH, the shape
functions failed.
In the attachment, I showed the information for the first active cell only.
It is easy to get information of other cells with the codes.

I think that there exists something wrong with the cell that are not a
cube, but a hexahedron when initializing the shape functions (fe_values).

I think it would be easier for you to help me with the attachment below.
Thank you very much.
Best regards.
Pham Ngoc Kien

Vào Th 2, 15 thg 4, 2019 vào lúc 21:48 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 4/11/19 9:17 PM, Phạm Ngọc Kiên wrote:
> > Testing for an edge whose global vertices located from (0,0,0) to
> (0,0,1) in
> > real coordinates.
> > With a cube I get the shape function vectors at the dof related to the
> edge,
> > for examples, (0,0,0), (0,0,-0.25), (0,0,-0.5), (0,0,-1), which are
> parallel
> > to the edge.
> > However, if I create a mesh with GMSH and check the shape function
> vectors,
> > some of them are: (-0.0868383, -0.428, 1), (-0.129933, -0.402772, 1),
> > (-0.052103 -0.2568 0.6),..,
> > which are not parallel to the edge.
> >
> > In my limited opinion, I will get the wrong solution as the shape
> functions
> > created by my code are not parallel to the edge in my model.
> > Could you please tell me why this happen and how to fix it?
>
> Good question. This does sound wrong. To figure out what is happening
> exactly,
> we typically construct small test cases that try to illustrate one
> particular
> issue. In your case, I would try to create a mesh with exactly one cell --
> either the one you have above, or maybe even simpler a cube that is just
> rotated. Then create a DoFHandler on it and output the values of the shape
> functions. This way, you really have only one thing that could go wrong,
> and
> it's easy to demonstrate what you see and how it differs from expectations.
>
> Do you think you could come up with a small program that does this?
>
> Best
>   WB
>
> --
> 
> 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.
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.
/* -
 *
 * Copyright (C) 2013 - 2018 by the deal.II authors
 *
 * This file is part of the deal.II library.
 *
 * The deal.II library is free software; you can use it, redistribute
 * it, and/or modify it under the terms of the GNU Lesser General
 * Public License as published by the Free Software Foundation; either
 * version 2.1 of the License, or (at your option) any later version.
 * The full text of the license can be found in the file LICENSE.md at
 * the top level directory of deal.II.
 *
 * -

 *
 * Author: Pham Ngoc Kien, Seoul National University, 2019
 */


#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 


#include 
#include 
#include 
#include 

#include 
#include 
#include 
#include 


#include 
#include 
#include 
#include 
#include 

#include 
#include 
#include 
#include 


#include 
#include 

#include 


#include 
#include 
#include 

#include 
#include 
#include 

using namespace dealii;

template
class EM {
public:
EM();

~EM();

void run();


private:
void grid_generation();

void setup_system();

void assemble_system();


void print_mesh_info(const Triangulation ,
 const std::string );


Triangulation triangulation;
MappingQ mappin

[deal.II] The shape function on physical cell in hexahedral is not parallel with the edge (FeNedelecSZ and FENedelec)

2019-04-11 Thread Phạm Ngọc Kiên

Hi colleagues,
When using FE_NedelecSZ in my code.
Testing for an edge whose global vertices located from (0,0,0) to (0,0,1) 
in real coordinates.
With a cube I get the shape function vectors at the dof related to the 
edge, for examples, (0,0,0), (0,0,-0.25), (0,0,-0.5), (0,0,-1), which are 
parallel to the edge.
However, if I create a mesh with GMSH and check the shape function vectors, 
some of them are: (-0.0868383, -0.428, 1), (-0.129933, -0.402772, 1),  
(-0.052103 -0.2568 0.6),.., 
which are not parallel to the edge.

In my limited opinion, I will get the wrong solution as the shape functions 
created by my code are not parallel to the edge in my model.
Could you please tell me why this happen and how to fix it?

I would like to thank you very much in advance.
Best regards,.

-- 
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.

Re: [deal.II] Why the FENedelec and FENedelecSZ 's shape functions change versus edge length?

2019-04-10 Thread Phạm Ngọc Kiên

Thank you very much.
Could you please recommend me some articles or books which I can find the
idea to compute the integral in real space for assembling system?
I think it would be really helpful for me.
Best regards,

Vào Th 5, 11 thg 4, 2019 vào lúc 12:49 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 4/10/19 7:59 PM, Phạm Ngọc Kiên wrote:
> >
> > In finite element method, we transform the integral from the physical to
> the
> > reference coordinate system.
> > Thus to compute  \int \varphi_i(x)  \varphi_j(x)  dx in a physical cell,
> I
> > will need:
> > \varphi_i(x_q), \varphi_j(x_q), and JxW in reference cell to compute the
> > integral by quadrature formula.
> > The term |det J(x_q)| is to change dx in physical coordinates into d(
> \hat x)
> > in reference coordinates so that we can compute the integral in
> reference cell.
> > However, in my code fe_values[vec[block_index_i]].value(i, q_point) is
> the
> > shape function in physical coordinate system.
>
> Correct. In deal.II, you do the integration in real space. You don't need
> to
> do the transformation back to the reference cell. (Of course, this happens
> somewhere internally to the library, but it is not something you need to
> worry
> about.)
>
>
> > If I use:
> > sum of quadrature point { fe_values[vec[block_index_i]].value(i,
> q_point) *
> > fe_values[vec[block_index_j]].value(j, q_point) *fe_values.JxW(q_point)
> }
> > It means that I am using the shape functions in physical coordinate
> system to
> > compute the integral instead of those in reference one.
>
> Yes. And that's correct.
>
> Best
>   WB
>
>
> --
> 
> 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.
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] Why the FENedelec and FENedelecSZ 's shape functions change versus edge length?

2019-04-10 Thread Phạm Ngọc Kiên

I am sorry about my unclear question.

In finite element method, we transform the integral from the physical to
the reference coordinate system.
Thus to compute  \int \varphi_i(x)  \varphi_j(x)  dx in a physical cell, I
will need:
\varphi_i(x_q), \varphi_j(x_q), and JxW in reference cell to compute the
integral by quadrature formula.
The term |det J(x_q)| is to change dx in physical coordinates into d( \hat
x) in reference coordinates so that we can compute the integral in
reference cell.
However, in my code fe_values[vec[block_index_i]].value(i, q_point) is the
shape function in physical coordinate system.
If I use:
sum of quadrature point { fe_values[vec[block_index_i]].value(i, q_point) *
fe_values[vec[block_index_j]].value(j, q_point)  *fe_values.JxW(q_point)  }
It means that I am using the shape functions in physical coordinate system
to compute the integral instead of those in reference one.
Thus, I think I am wrong here, isn't it?

As my shape function is a vector function, its values in reference
coordinates are different from those in physical coordinates.
However, the FE_NedelecSZ and FE_Nedelec classes do not provide
fe.shape_value() and fe.shape_value_component() because the shape function
is not primitive.
Therefore, I am trying to get the shape function in reference cell to put
in my code.

I hope that this question is clear enough for you to help me.
Thank you very much.

Best regards,
Pham Ngoc Kien

Vào Th 4, 10 thg 4, 2019 vào lúc 12:56 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 4/9/19 8:51 PM, Phạm Ngọc Kiên wrote:
> > I am a little bit confusing with the integral over all cell :
> > *\int \varphi_i(x) *\varphi_j(x) * dx *on physical cell is approximated
> > by computing *\sum { \varphi_hat_i(x_hat) *\varphi_hat_i(x_hat)
> *JxW(x_hat) }
> > *on reference cell.
> > With FEValues we get *\varphi_i(x) *, if we compute the above integral
> as:
> > *\sum { \*varphi_i(x)* *\*varphi_j(x)* *JxW(x_hat) }*
> > then the value of the integral changes because shape functions defined
> on
> > reference cell are different from those on real cell.
> >
> > Is there any reason for the integral calculation*\sum { \*varphi_i(x)*
> > *\*varphi_j(x)* *JxW(x_hat) }?*
>
> I have to admit that I don't understand the question. The integral you are
> computing is
>
>\int \varphi_i(x)  \varphi_j(x)  dx
>
> which is approximated using quadrature via
>
>\sum_q  \varphi_i(x_q) \varphi_j(x_q) w_q
>
> where w_q = |det J(x_q)| \hat w_q = JxW. So this approximation happens in
> real
> space, not on the reference cell. How the individual components are
> computed
> (via the reference cell) is not important for the approximation.
>
> Best
>   WB
>
> --
> 
> 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.
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] Why the FENedelec and FENedelecSZ 's shape functions change versus edge length?

2019-04-09 Thread Phạm Ngọc Kiên

I am a little bit confusing with the integral over all cell :
*\int \varphi_i(x)  \varphi_j(x)  dx *on physical cell is approximated
by computing *\sum { \varphi_hat_i(x_hat)  \varphi_hat_i(x_hat) JxW(x_hat)
} *on reference cell.
With FEValues we get *\varphi_i(x) *, if we compute the above integral as:
*\sum { \varphi_i(x) \varphi_j(x) JxW(x_hat) }*
then the value of the integral changes because shape functions defined on
reference cell are different from those on real cell.

Is there any reason for the integral calculation* \sum { \varphi_i(x)
\varphi_j(x) JxW(x_hat) }?*
Thank you very much for your quick answers.

Best regards.


Vào Th 4, 10 thg 4, 2019 vào lúc 00:15 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 4/9/19 1:27 AM, Phạm Ngọc Kiên wrote:
> > I looked through the codes for Fe_Nedelec_SZ in the library and found
> that the
> > fe_values was transformed from reference cell to physical cell by the
> function
> > like:
> > mapping.transform(make_array_view(fe_data.shape_values[dof]),
> >  mapping_covariant,
> >  mapping_internal,
> >  make_array_view(transformed_shape_values));
> > Does this mean that fe_values.value(i,q_point) and
> fe_values.curl(i,q_point)
> > are stored on physical cell?
>
> Yes, yes, of course -- that's the point of FEValues :-)
>
> 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.
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] Why the FENedelec and FENedelecSZ 's shape functions change versus edge length?

2019-04-09 Thread Phạm Ngọc Kiên

Dear all,
I looked through the codes for Fe_Nedelec_SZ in the library and found that
the fe_values was transformed from reference cell to physical cell by the
function like:
mapping.transform(make_array_view(fe_data.shape_values[dof]),
mapping_covariant,
mapping_internal,
make_array_view(transformed_shape_values));
Does this mean that fe_values.value(i,q_point) and
fe_values.curl(i,q_point) are stored on physical cell?

In mapping, mapping of integrals
<https://www.dealii.org/current/doxygen/deal.II/classMapping.html#a3ce74638006ce574195c55a69aa48d79>,
the shape function is defined on reference cell to approximate the integral
by quadrature.
How can I use fe_values.value(i,q_point) and fe_values.curl(i,q_point) to
compute the integral when they are defined on physical cell?
Is there any other way to compute the integral in this case?

I would like to thank you very much in advance.
Best regards,
Pham Ngoc Kien


Vào Th 5, 28 thg 3, 2019 vào lúc 16:30 Phạm Ngọc Kiên <
ngockien.lea...@gmail.com> đã viết:

> Thank you very much
> I will dig deeper from your suggestion.
> Best regards,
> Pham Ngoc Kien
>
> Vào Th 5, 28 thg 3, 2019 vào lúc 12:25 Wolfgang Bangerth <
> bange...@colostate.edu> đã viết:
>
>> On 3/27/19 6:56 PM, Phạm Ngọc Kiên wrote:
>> > Yes, I know that the source should be integral over the edge of a cell,
>> not a
>> > single point on the cell.
>> > And I checked that if I take the integral over an edge of shape
>> function I
>> > will get 1 * inverse of edge length.
>> > However, we actually do compute the source term over a cell when
>> assembling
>> > the right hand side:
>> > F_i  =  \int_\Omega  \phi_i(x) J(x)  dx
>> > To compute it, we implement the Quadrature formula with its weights on
>> unit
>> > cell rather than the weights on and edge.
>> > In my limited knowledge, I think that J(x) should be non-zero if x is
>> on my
>> > source line, and be zero every where else in the cell.
>>
>> Well, but that's precisely a delta function then. J is a current
>> *density*, so
>> if you have a current density that only lives on an edge with a finite
>> length
>> but an infinitely small diameter, then the total current that is flowing
>> is
>> zero. If you say that J(x) is zero everywhere but on that edge, then no
>> current is flowing. Your model is just flawed in that case.
>>
>> If you want to restrict a *nonzero* current to just an edge, then you
>> need an
>> infinite current density, and in that case you would have
>>
>>F_i = \int_Omega phi_i(x) J(x)  // J is current *density*
>>= \int_edge phi_i(x) j(x)   // j is *current*
>>
>> where the second integral only extends over a one-dimensional edge.
>>
>>
>> > Should I set J(x) constant when computing F_i if dof i of the cell
>> relates to
>> > the source?
>>
>> You're asking the wrong person. If you are trying to model a physical
>> situation, then what J(x) is is a question of how that physical situation
>> looks like and is best described. Nobody but you can answer this -- it's
>> not a
>> mathematical question, but a modeling question.
>>
>> 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.
>> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] Why the FENedelec and FENedelecSZ 's shape functions change versus edge length?

2019-03-28 Thread Phạm Ngọc Kiên

Thank you very much
I will dig deeper from your suggestion.
Best regards,
Pham Ngoc Kien

Vào Th 5, 28 thg 3, 2019 vào lúc 12:25 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 3/27/19 6:56 PM, Phạm Ngọc Kiên wrote:
> > Yes, I know that the source should be integral over the edge of a cell,
> not a
> > single point on the cell.
> > And I checked that if I take the integral over an edge of shape function
> I
> > will get 1 * inverse of edge length.
> > However, we actually do compute the source term over a cell when
> assembling
> > the right hand side:
> > F_i  =  \int_\Omega  \phi_i(x) J(x)  dx
> > To compute it, we implement the Quadrature formula with its weights on
> unit
> > cell rather than the weights on and edge.
> > In my limited knowledge, I think that J(x) should be non-zero if x is on
> my
> > source line, and be zero every where else in the cell.
>
> Well, but that's precisely a delta function then. J is a current
> *density*, so
> if you have a current density that only lives on an edge with a finite
> length
> but an infinitely small diameter, then the total current that is flowing
> is
> zero. If you say that J(x) is zero everywhere but on that edge, then no
> current is flowing. Your model is just flawed in that case.
>
> If you want to restrict a *nonzero* current to just an edge, then you need
> an
> infinite current density, and in that case you would have
>
>F_i = \int_Omega phi_i(x) J(x)  // J is current *density*
>= \int_edge phi_i(x) j(x)   // j is *current*
>
> where the second integral only extends over a one-dimensional edge.
>
>
> > Should I set J(x) constant when computing F_i if dof i of the cell
> relates to
> > the source?
>
> You're asking the wrong person. If you are trying to model a physical
> situation, then what J(x) is is a question of how that physical situation
> looks like and is best described. Nobody but you can answer this -- it's
> not a
> mathematical question, but a modeling question.
>
> 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.
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] Why the FENedelec and FENedelecSZ 's shape functions change versus edge length?

2019-03-27 Thread Phạm Ngọc Kiên

Dear Prof. Bangerth,
Yes, I know that the source should be integral over the edge of a cell, not
a single point on the cell.
And I checked that if I take the integral over an edge of shape function I
will get 1 * inverse of edge length.
However, we actually do compute the source term over a cell when assembling
the right hand side:
F_i  =  \int_\Omega  \phi_i(x) J(x)  dx
To compute it, we implement the Quadrature formula with its weights on unit
cell rather than the weights on and edge.
In my limited knowledge, I think that J(x) should be non-zero if x is on my
source line, and be zero every where else in the cell.
Should I set J(x) constant when computing F_i if dof i of the cell relates
to the source?

Now I am considering an edge rather than a cell in the whole domain is the
source.
Could you please give me an idea?

I would like to thank you very much for your answers.

Best regards,

Vào Th 4, 27 thg 3, 2019 vào lúc 22:25 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 3/27/19 2:52 AM, Phạm Ngọc Kiên wrote:
> >
> > My question for assembling the system is
> > when I set the cell right hand side is none zero at a given dofs (source
> > dofs), should I assemble the cell matrix at these dofs none zeros?
> > Should both cell matrix (i,j) and cell right hand side (i) be zeros for
> the
> > dofs i , j  which are not at the source dof?
>
> I don't think we quite understood what you are doing. Let's say your right
> hand side function is J(x), then the right hand side vector should be
>
>F_i  =  \int_\Omega  \phi_i(x) J(x)  dx
>
> I *think* what you are describing is that your J(x) is actually a delta
> function at a single point, which you call the "source dof". If that is
> the
> case, say that J(x) is a delta function located at x0, then you need to
> compute
>
>F_i  =  \phi_i(x0)
>
> which may or may not be equal to one -- in any case, making the assumption
> that a shape function has a particular value at a particular point is a
> bad
> idea, and the right approach is to evaluate it at that point to get its
> value.
>
> It is worth pointing out that while the description of source located at
> individual points is often used, it is not a physically correct
> description.
> No sources are point sources. A physically correct description would use
> distributed sources (which leads to integrals that are not approximated by
> point values), and as you are discovering, using the "simplification" of
> describing sources as point sources actually makes all sorts of things
> more
> complicated in finite element codes.
>
> Does this make sense?
>
> 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.
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] Why the FENedelec and FENedelecSZ 's shape functions change versus edge length?

2019-03-27 Thread Phạm Ngọc Kiên

Dear Prof. Bangerth,
In my limited knowledge, I think that my solution wrong come from the
assembling the matrix and right hand side.
Because as you mentioned, the solution remains the same no matter how I
change the mesh.

My question for assembling the system is
when I set the cell right hand side is none zero at a given dofs (source
dofs), should I assemble the cell matrix at these dofs none zeros?
Should both cell matrix (i,j) and cell right hand side (i) be zeros for the
dofs i , j  which are not at the source dof?

I saw from testing my codes that the cell right hand side also changed when
the volume of the cell changed although I remained the length of 1 edge.
Is this because the integral over all cell also be scaled by a certain
factor?
Thank you very much.
 Best regards,
Pham Ngoc Kien

Vào Th 3, 26 thg 3, 2019 vào lúc 07:48 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 3/22/19 5:39 PM, Phạm Ngọc Kiên wrote:
> > 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?
>
> Well, the "scaling" is obtained by computing the solution as
>
>u_h(x) = \sum_j  U_j  \varphi_j(x)
>
> The values of U_j by itself are meaningless without the shape functions
> \varphi_j. If someone decided to multiply every shape function by ten,
> but then also divided all of the U_j by ten, then the solution u_h(x)
> remains the same.
>
> 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.
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] Why the FENedelec and FENedelecSZ 's shape functions change versus edge length?

2019-03-22 Thread Phạm Ngọc Kiên

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
void OutputTools::point_solution(const Mapping ,
 const DoFHandler _handler,
 const Vector ,
 const Point ,
 Vector _solution) const {

const std::pair::active_cell_iterator, Point>
cell_point =
GridTools::find_active_cell_around_point(mapping, dof_handler, point);

const double delta = 1e-12;
Assert(GeometryInfo::distance_to_unit_cell(cell_point.second) < delta,
   ExcInternalError());

const Quadrature
quadrature(GeometryInfo::project_to_unit_cell(cell_point.second));
FEValues 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> temp_solution_re(quadrature.size());
std::vector> 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 values (dim*2);
VectorTools::point_value(dof_handler,solution,point,values);
std::cout<<"\n";
std::cout<<"cell: \t"< đã 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 dealii+unsubscr...@googlegroups.com.
> 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"

Re: [deal.II] Why the FENedelec and FENedelecSZ 's shape functions change versus edge length?

2019-03-21 Thread Phạm Ngọc Kiên

Dear Prof.Wolfgang,
I checked again my codes and looked at the shape functions on the edge.
When I looked at one edge of a cell with its degree of freedom, for
example, with the edge length is 1, I saw the shape functions at quadrature
points on the edge is (0,0,1).
Thus, if I calculate the integral over the edge by hand and  got the result
is 1 as what you said.

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.

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?

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.

I think there should be some where I was wrong when assembling the matrix
and the right hand side.

I would like to thank you very much for your great support.
Best regards,
Pham Ngoc Kien



Vào Th 5, 21 thg 3, 2019 vào lúc 22:54 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

>
> Pham,
>
> > When  testing my codes with FENedelec and FENedelecSZ elements, I saw
> that the
> > shape functions change versus the mesh size.
> > I think that for each edge of a cell, the shape functions for the degree
> of
> > freedom related to that edge are scaled with the inverse edge length.
> > For example, assuming I have a cell that is a cube, whose edge length is
> 1,
> > the shape function at the quadrature point located on the edge is {0, 0,
> 1}.
> > For the smaller cell whose edge length is 0.5, the shape function values
> is
> > {0, 0, 2}.
>
> Yes. This can be because of two things:
> * Shape functions are transformed by the Piola or a related transform that
>multiplies shape functions by a measure related to the size of the cell.
> * Shape functions can be defined not as point values at specific points,
>but as *integrals* over an edge or a face. So if the degree of freedom
>is that the integral of a shape function (or its tangential component)
>over an edge needs to be 1.0, then the *height* of the function needs to
>increase by a factor of 2 if the size of the cell is reduced by a
> factor of
>2.
>
>
> > In my limited opinion, the shape function is defined in unit cell over
> > [0,1]^3, and the shape function should not change when the real cell
> change
> > its size.
>
> This is not necessarily so. It depends on the finite element and the
> corresponding transform, see above.
>
>
> > This problem leads to the change of my solution when using the finer
> mesh
> > compared to that of coarser one.
>
> What is the problem you are encountering?
>
> If you refine the mesh and the height of the shape functions increase by a
> factor, then the entries of the solution vector will just go down by a
> factor
> of 2. The product of shape function times solution coefficients will then
> remain the same.
>
> 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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Why the FENedelec and FENedelecSZ 's shape functions change versus edge length?

2019-03-19 Thread Phạm Ngọc Kiên

Hi colleagues,
When  testing my codes with FENedelec and FENedelecSZ elements, I saw that 
the shape functions change versus the mesh size.
I think that for each edge of a cell, the shape functions for the degree of 
freedom related to that edge are scaled with the inverse edge length.
For example, assuming I have a cell that is a cube, whose edge length is 1, 
the shape function at the quadrature point located on the edge is {0, 0, 1}.
For the smaller cell whose edge length is 0.5, the shape function values is 
{0, 0, 2}. 

In my limited opinion, the shape function is defined in unit cell over 
[0,1]^3, and the shape function should not change when the real cell change 
its size.
This problem leads to the change of my solution when using the finer mesh 
compared to that of coarser one. 
Because both my cell matrix and cell rhs change versus the edge length. 

For FeNedelecSZ, I tested with a cell limited in the region [0,1]^3, which 
means that the edge length is 1 everywhere.
When printing out the integral over the cell for one degree of freedom by 
   
\int_Omega  \varphi_i(x)  Js(x)  dx 
For simplicity, I set Omega =1 and Js(x) = {1,1,1} for testing the integral 
over the cell at degree of freedom i_th when source magnitude is 1.
The result for this is 0.25 rather than 1.
Is this means that for a source of magnitude 1., I actually can model it as 
a source edge at degree of freedom i_th with the magnitude 0.25?

I think there should be some idea that we change the shape function values 
like this.
Could you please give me some advice?

Thank you very much for your great support.
Best regards,

-- 
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.

Re: [deal.II] what is the different between VectorTools::point_value and fe_values.get_function_values()?

2019-03-18 Thread Phạm Ngọc Kiên

Thank you very much for your guidance.
I think that with the H curl conforming element like FeNedelec or
FeNedelecSZ, my solution vector has continuous tangential component, isn't
it?
Thus, Is there some where I can see my vector solution is discontinuous?
Can I check it by visualization?

I think my below question is not related to the topic, but may I ask some
thing about the shape function of FeNedelec or FeNedelecSZ?
My testing codes showed that with those elements, the shape functions
changed the values when the mesh grid changed.
For example, assuming I have a cell that is a cube, whose edge length is 1,
the shape function at the edge is {0, 0, 1}.
For the smaller cell whose edge length is 0.5, the shape function at the
edge is {0, 0, 2}.
I saw it was double every time I reduce the edge length by the factor of 2.
In my limited opinion, the shape function is defined in unit cell over
[0,1]^3, and the shape function should not change when the real cell change
its size.
This problem leads to my solution change when refining the mesh.
I think there should be some idea that we change the shape function values
like this.
Could you please give me some advice?

Best regards,
Pham Ngoc Kien

Vào Th 3, 19 thg 3, 2019 vào lúc 02:22 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 3/18/19 4:14 AM, Daniel Arndt wrote:
> >
> > VectorTools::point_value() does more or less the same you are doing
> above
> > manually, i.e. calling GridTools::find_active_cell_around_point() and
> then
> > initializing a FEValues object for evaluating the given finite element
> vector.
>
> Specifically, you can take a look here:
>
> https://github.com/dealii/dealii/blob/master/include/deal.II/numerics/vector_tools.templates.h#L8633
>
> 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.
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] what is the different between VectorTools::point_value and fe_values.get_function_values()?

2019-03-14 Thread Phạm Ngọc Kiên

I highly appreciate all of your answers. It becomes clearly for me to look
through the library.
Let me show all of my function and my ideas below.
I set the point for evaluating its solution to be on the vertex of a cell,
and using the following function to get the solution

template
void OutputTools::point_solution(const Mapping ,
//input mapping
  const DoFHandler
_handler,   //input dof_handler
  const Vector ,
 //input solution vector
  const Point ,
 //input point in real coordinates to get the solution
  Vector _solution)
const {  //output

const std::pair::active_cell_iterator, Point>
cell_point =
GridTools::find_active_cell_around_point(mapping, dof_handler, point);
 //find active cell around input point

const double delta = 1e-20;
Assert(GeometryInfo::distance_to_unit_cell(cell_point.second) < delta,
   ExcInternalError());

const Quadrature
quadrature(GeometryInfo::project_to_unit_cell(cell_point.second));
//get quadrature by project the point to unit cell
FEValues 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> temp_solution_re(quadrature.size());
std::vector> 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];
}

   //using point_value() for comparing the result

   Vector values (6);
   VectorTools::point_value(dof_handler,solution,point,values);

//print out all the results for testing the code

   std::cout<<"\n";
   std::cout<<"cell: \t"<::active_cell_iterator, Point>
cell_point =
GridTools::find_active_cell_around_point(mapping, dof_handler, point);

I personally think that as I set it to be const, so the return cell
(cell_point.first) is only one cell, rather than a group of cells.
Could you please tell me about this?

2. My second idea is with the codes:
const Quadrature
quadrature(GeometryInfo::project_to_unit_cell(cell_point.second));
FEValues fe_values(mapping, dof_handler.get_fe(), quadrature,
update_values | update_gradients |
update_quadrature_points);
fe_values.reinit(cell_point.first);

I have only one quadrature defined on the cell, and the printing
result show me the same thing.
Can I use this way to get the solution, or it would be more efficient
when using codes like:
QGaussLobatto quadrature_formula(2);
FEValues fe_values(fe, quadrature_formula,
update_values | update_gradients |
update_quadrature_points | update_JxW_values);
fe_values.reinit(cell_point.first);

And then finding the quadrature point coincide with the evaluating
point to get the solution?

3. My third thought is in my function I have only 1 quadrature point
so that temp_solution_re[0] and temp_solution_im[0] are the variables
containing the solution.

The printing results showed that the solution from point_value() is
identical with those from temp_solution_re[0] and temp_solution_im[0].

Do temp_solution_re[1], temp_solution_re[2],
temp_solution_im[1],temp_solution_im[2] make any sense for my
solution? Because I saw the code can run with these variables.

Thank you very much!

Vào Th 6, 15 thg 3, 2019 vào lúc 00:29 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 3/14/19 1:28 AM, Phạm Ngọc Kiên wrote:
> >
> > So for the first function
> > VectorTools::point_value(dof_handler,solution,point,values);
> > I can specify the observed point in real coordiante, solution vector and
> > dof_handler to get the values, right?
>
> Correct.
>
>
> > For the second one, assuming that I am using QGauss for solving the
> > system, does it mean that if the point is on the edge of the cell then I
> > cannot get the solution at that point?
>
> Yes. It can be any quadrature formula defined on the reference cell. It
> doesn't matter where the points are located.
>
> 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/foru

Re: [deal.II] what is the different between VectorTools::point_value and fe_values.get_function_values()?

2019-03-14 Thread Phạm Ngọc Kiên

Thank you very much for your answer.

So for the first function  
VectorTools::point_value(dof_handler,solution,point,values); 

I can specify the observed point in real coordiante, solution vector and 
dof_handler to get the values, right?

For the second one, assuming that I am using QGauss for solving the system, 
does it mean that if the point is on the edge of the cell then I cannot get 
the solution at that point?
I am using the codes like this:

const std::pair::active_cell_iterator, Point>
cell_point = GridTools::find_active_cell_around_point(mapping, 
dof_handler, point);

const double delta = 1e-20;
Assert(GeometryInfo::distance_to_unit_cell(cell_point.second) < delta,
   ExcInternalError());

const Quadrature 
quadrature(GeometryInfo::project_to_unit_cell(cell_point.second));
FEValues 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> temp_e_field_re(1);
std::vector> temp_e_field_im(1);

//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);


With this function, I see that it returns 9 components: 3 components for 
each parameters: temp_solution_re[0], temp_solution_re[1], temp_solution_re[
2] 
when I print them out by:

std::cout<<"result with get_function_values real: \t"<< temp_solution_re[0] << 
", "<
> On 3/13/19 8:21 PM, Phạm Ngọc Kiên wrote: 
> > I am testing my codes for output the solution at a point in my numerical 
> model. 
> > I saw in the library 2 ways to do this task. The first one is to use 
> > VectorTools::point_value() function and the second one is 
> > fe_values.get_function_values(). 
>
> These functions are really used in very different contexts. The first of 
> these 
> is when you want to evaluate a solution vector at some arbitrary point in 
> the 
> domain. As a consequence, you have to specify this point when you call the 
> function. 
>
> The second function is used when you want to evaluate the function at very 
> specific points, namely the quadrature points of the current cell. You can 
> of 
> course only use this if the point at which you want to evaluate the 
> solution 
> happens to be a quadrature point of the current cell, and not just any 
> point 
> anywhere. In return the second function is then vastly more efficient than 
> the 
> first. 
>
> Does this help? 
>
> 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.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] what is the different between VectorTools::point_value and fe_values.get_function_values()?

2019-03-14 Thread Phạm Ngọc Kiên

Thank you very much for your answer.
So for the first function
VectorTools::point_value(dof_handler,solution,point,values);

I can specify the observed point in real coordiante, solution vector and
dof_handler to get the values right?

For the second one, assuming that I am using QGauss for solving the system,
can I use QGaussLobatto to get the quadrature points on the edge of a cell
and get the values from this?
With the second function, I see that it returns 9 components: 3 for each
parameters: temp_solution_re[0], temp_solution_re[1], temp_solution_re[3]
when I call fe_values[vector_re].get_function_values(solution,
temp_solution_re);
I really want to know which one is the solution for my vector in x,y,z
direction?

Best regards,
Pham Ngoc Kien

Vào Th 5, 14 thg 3, 2019 vào lúc 14:08 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 3/13/19 8:21 PM, Phạm Ngọc Kiên wrote:
> > I am testing my codes for output the solution at a point in my numerical
> model.
> > I saw in the library 2 ways to do this task. The first one is to use
> > VectorTools::point_value() function and the second one is
> > fe_values.get_function_values().
>
> These functions are really used in very different contexts. The first of
> these
> is when you want to evaluate a solution vector at some arbitrary point in
> the
> domain. As a consequence, you have to specify this point when you call the
> function.
>
> The second function is used when you want to evaluate the function at very
> specific points, namely the quadrature points of the current cell. You can
> of
> course only use this if the point at which you want to evaluate the
> solution
> happens to be a quadrature point of the current cell, and not just any
> point
> anywhere. In return the second function is then vastly more efficient than
> the
> first.
>
> Does this help?
>
> 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.
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[deal.II] what is the different between VectorTools::point_value and fe_values.get_function_values()?

2019-03-13 Thread Phạm Ngọc Kiên

Hi colleagues,
I am testing my codes for output the solution at a point in my numerical 
model.
I saw in the library 2 ways to do this task. The first one is to use 
VectorTools::point_value() function and the second one is 
fe_values.get_function_values().
For the latter function, I am using the code like this:

std::vector> temp_solution_re(1); 
std::vector> temp_solution_im(1);  //as my solution have both 
real and imaginary parts
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) {
solution(j) = temp_solution_re[0][j];
solution(j+dim) = temp_solution_im[0][j];
}
Which returns the resutls are a vector of 3 component in x,y,z of my 
solution,and is the same as using
Vector values (dim+dim);
VectorTools::point_value(dof_handler,solution,point,values);

However, I checked that each tensor temp_solution_re[0], temp_solution_re[1], 
temp_solution_re[3] has 3 component in x,y,z direction.
Could you please tell me the meaning of these tensors when we get the solution 
by using fe_values.get_function_values() ?
Thank you very much.

Best regards,
Pham Ngoc Kien

-- 
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.

Re: [deal.II] How to get the coordinates of a given degree of freedom on an edge?

2019-02-22 Thread Phạm Ngọc Kiên

Yes.
Thank you very much.
I have just started coding with C++ in several month and this is really
helpful for me.
Best,
Kien

Vào Th 6, 22 thg 2, 2019 vào lúc 21:15 Jean-Paul Pelteret <
jppelte...@gmail.com> đã viết:

> Dear Kien,
>
> The problem that you’re facing (in the example that you gave us) is not
> one to do with deal.II itself, but rather C++ and how you’ve structured
> your code. When I build your code in debug mode, I see this warning:
>
> $ make
> Scanning dependencies of target problem
> [ 50%] Building CXX object CMakeFiles/problem.dir/problem.cc.o
> /problem.cc:112:31: warning: field 'order' is uninitialized when
> used here [-Wuninitialized]
>  fe(FE_NedelecSZ(order), 1,   //(order), multiplicity
>
> and when I run it, it actually hangs on my machine. What’s happening here
> relates to where you have this line:
>
> const unsigned int order = 3;
>
> It is in the body of the EM class, but most critically its below the
> declarations for
>
> DoFHandler dof_handler;
> FESystem fe;
>
> This means that the “order” variable is initialised AFTER the
> “dof_handler", and “fe”, and so during the call to the constructor it
> uninitialised when used to initialise the “dof_handler” and “fe”. So why
> this hangs on my machine is likely because “order” is random number, and
> more specifically something very large when I tested it. Its luck that it
> the section of memory assigned to this variable was either all zeros (with
> the interpretation of the uninitialised “order” variable being equal to
> zero, which corresponds to the lowest order Nedelec element that is
> permitted) or with some random bits that were interpreted as a small number.
>
> The quick fix to the issue is to move the definition of “order” above that
> of “dof_handler”. Then the warning goes away and everything works as
> expected. If I change the constructor definition slightly to
>
> template
> EM::EM()
> :dof_handler(triangulation),
>  fe(FE_NedelecSZ(order), 1,   //(order), multiplicity
> FE_NedelecSZ(order), 1)
> {
>   std::cout << "order: " << order << "  fe.dofs_per_cell: " <<
> fe.dofs_per_cell << std::endl;
> }
>
> then I get the following output for two different orders.:
> order: 2  fe.dofs_per_cell: 288
> order: 3  fe.dofs_per_cell: 600
>
> A better fix might be to pass the order into the constructor, or to use a
> parameter handler to choose the order at run-time.
>
> I hope that this helps you understand where things were going wrong.
>
> Best,
> Jean-Paul
>
> On 22 Feb 2019, at 02:37, Phạm Ngọc Kiên 
> wrote:
>
> Hi,
> The attachment is my small program that I created.
> My fe system have 24 dofs for both real and imaginary of 12 edges of a
> cell.
> No matter how I change the parameter "order" the number of dofs is
> unchanged.
> I think I have problem when constructing the object.
> I would like to thank you very much.
>  Best,
> Kien
>
> Vào Th 5, 21 thg 2, 2019 vào lúc 23:38 Wolfgang Bangerth <
> bange...@colostate.edu> đã viết:
>
>> On 2/20/19 7:01 PM, Phạm Ngọc Kiên wrote:
>> > I tested in my codes using
>> >
>> > template
>> > CSEM::CSEM()
>> >  ://mapping(1),
>> > dof_handler(triangulation),
>> >  fe(FE_NedelecSZ(order),1,//(order), multiplicity
>> > FE_NedelecSZ(order),1) {}
>> >
>> > When I changed the parameter order, I saw only the number of dofs was
>> only of lowest order dofs no matter how big the order was
>> >
>> > I think that there should be a way to get higher order dofs, but I
>> cannot find it now.
>>
>> Can you create a small program that demonstrates this issue?
>>
>> 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.
>> 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 opti

Re: [deal.II] Re: How to get the coordinates of a given degree of freedom on an edge?

2019-02-21 Thread Phạm Ngọc Kiên

Hi,
The attachment is my small program that I created.
My fe system have 24 dofs for both real and imaginary of 12 edges of a cell.
No matter how I change the parameter "order" the number of dofs is
unchanged.
I think I have problem when constructing the object.
I would like to thank you very much.
 Best,
Kien

Vào Th 5, 21 thg 2, 2019 vào lúc 23:38 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 2/20/19 7:01 PM, Phạm Ngọc Kiên wrote:
> > I tested in my codes using
> >
> > template
> > CSEM::CSEM()
> >  ://mapping(1),
> > dof_handler(triangulation),
> >  fe(FE_NedelecSZ(order),1,//(order), multiplicity
> > FE_NedelecSZ(order),1) {}
> >
> > When I changed the parameter order, I saw only the number of dofs was
> only of lowest order dofs no matter how big the order was
> >
> > I think that there should be a way to get higher order dofs, but I
> cannot find it now.
>
> Can you create a small program that demonstrates this issue?
>
> 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.
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.
/* -
 *
 * Copyright (C) 2013 - 2018 by the deal.II authors
 *
 * This file is part of the deal.II library.
 *
 * The deal.II library is free software; you can use it, redistribute
 * it, and/or modify it under the terms of the GNU Lesser General
 * Public License as published by the Free Software Foundation; either
 * version 2.1 of the License, or (at your option) any later version.
 * The full text of the license can be found in the file LICENSE.md at
 * the top level directory of deal.II.
 *
 * -


 */

#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 

#include 
#include 
#include 

#include 
#include 
#include 
#include 

#include 
#include 
#include 


#include 
#include 
#include 
#include 
#include 

#include 
#include 
#include 
#include 


#include 
#include 

#include 


#include 
#include 
#include 

#include 
#include 
#include 

using namespace dealii;

template
class EM {
public:
EM();


~EM();

void run();


private:
void grid_generation();

void setup_system();

void assemble_system();

void solve();



Triangulation triangulation;

DoFHandler dof_handler;
FESystem fe;

AffineConstraints constraints;

SparsityPattern sparsity_pattern;
SparseMatrix system_matrix;

Vector solution;
Vector system_rhs;

const unsigned int order = 3; //order of edge-based element

};

template
EM::EM()
:dof_handler(triangulation),
 fe(FE_NedelecSZ(order), 1,   //(order), multiplicity
FE_NedelecSZ(order), 1) {}

template
EM::~EM() {
dof_handler.clear();
}

template
void EM::grid_generation() {
GridGenerator::hyper_cube(triangulation, -4, 4);
}

template
void EM::setup_system() {
dof_handler.distribute_dofs(fe);

solution.reinit(dof_handler.n_dofs());
system_rhs.reinit(dof_handler.n_dofs());

constraints.clear();
DoFTools::make_hanging_node_constraints(dof_handler, constraints);

VectorTools::project_boundary_values_curl_conforming_l2(dof_handler,
0,
ConstantFunction(0., fe.n_components()),
1,  
constraints);
VectorTools::project_boundary_values_curl_conforming_l2(dof_handler,
dim,
ConstantFunction(

Re: [deal.II] Re: How to get the coordinates of a given degree of freedom on an edge?

2019-02-20 Thread Phạm Ngọc Kiên

Dear Daniel,
I think that I can address my problem with the source now after looking
through the way we create dofs in FE_NedelecSZ.
I tested in my codes using

template
CSEM::CSEM()
://mapping(1),
dof_handler(triangulation),
fe(FE_NedelecSZ(order), 1,   //(order), multiplicity
   FE_NedelecSZ(order), 1) {}

When I changed the parameter order, I saw only the number of dofs was
only of lowest order dofs no matter how big the order was

I think that there should be a way to get higher order dofs, but I
cannot find it now.

Thank you very much.


Vào Th 4, 20 thg 2, 2019 vào lúc 18:54 Daniel Arndt <
daniel.ar...@iwr.uni-heidelberg.de> đã viết:

> Kien,
>
> I am not quite sure I understand how the source term you want to use looks
> like.
>
> If you want to describe a point source, you should gave a look at
> https://www.dealii.org/developer/doxygen/deal.II/namespaceVectorTools.html#ac4ffc41e72e846029d04e7034de60d09
> .
> Be aware of the last paragraphs in its documentation, though.
>
> On the other hand, if you have a usual right hand side Js(x), you need to
> follow the approach Wolfgang described and that is used in almost all
> tutorial progams,
> e.g.
> https://www.dealii.org/developer/doxygen/deal.II/step_4.html#Step4assemble_system
> .
>
> How does Js(x) look like?
>
> I also tried to use FE_NedelecSZ(order) with different order.
>> And I saw that no matter how I change the order value, the functions
>> fe.n_dofs_per_quad() and fe.n_dofs_per_hex() returned 0. The number of dofs
>> was unchanged and was equal to fe.n_dofs_per_line().
>> Does this mean that there is no face and cell dof?
>>
> According to the implementation there should be
> - 0 dofs per vertex
> - degree+1 dofs per line
> - 2*degree*(degree+1) dofs per quad
> - 3*degree*degree*(degree+1) dofs per hex
>
>
>> Do I need to use hp::DoFHandler to get the higher order dofs for
>> calculation?
>>
> No, you would use this if you want to use a different ansatz space, e.g.
> the same type of element with a different polynomial degree, for each cell.
>
>  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.
>

-- 
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.

Re: [deal.II] Re: How to get the coordinates of a given degree of freedom on an edge?

2019-02-19 Thread Phạm Ngọc Kiên

Yes.
Let me go from my basic idea when setting up my source term.
First, I create a computational domain such that a source can be placed
exactly on an edge. It is always possible to do with unstructured meshes.
Then I go through each element and check if its edge overlap the source. If
it does then I will assign a source moment to the edge.
My problem now is to find out if any edge of an element is overlapped by
the my source.

I thought about finding the position of the degree of freedom in a cell (in
real co-ordinate) and check if its real position is my source position.
Thus, I sent a question about this.

I have checked my codes and saw that \varphi_i(x) were a vector at
quadrature point x for the i_th degree of freedom.
It determines the direction of my vector field on the i_th degree of
freedom in reference cell.
So, I set Js(x) is a vector of the same size to compute the right hand side.

I think there is no problem for the FE_NedelecSZ element. My code still run
but it do not give me an exact solution compared with analytic one because
it is not true if I set the same Js(x) for every dof in a cell containing
the source point.

I also tried to use FE_NedelecSZ(order) with different order.
And I saw that no matter how I change the order value, the functions
fe.n_dofs_per_quad() and fe.n_dofs_per_hex() returned 0. The number of dofs
was unchanged and was equal to fe.n_dofs_per_line().
Does this mean that there is no face and cell dof?
Do I need to use hp::DoFHandler to get the higher order dofs for
calculation?

Thank you very much.
Best,
Kien

Vào Th 4, 20 thg 2, 2019 vào lúc 07:44 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 2/19/19 12:03 AM, Phạm Ngọc Kiên wrote:
> > I currently try to solve the electromagnetics problem, followed by the
> > research of Ross Kynch.
> > My goal is to compute the right hand side for the problem:
> > curl curl E + i*omega*mu*sigma*E = - Js
> > I only want to set the Js term here to be none zero, for example to be
> > 1, at a given degree of freedom (i.e. at source position), and  0
> > everywhere else.
> > Thus, I need to write the codes to query if a degree of freedom is at
> > source position.
> > As you mentioned, now I think that I need to check if the source point
> > is on an edge when assembling the right hand side vector.
> > Could you please tell me help me about this?
>
> The term you are trying to add (Js) is a right hand side term to the
> equation. So you need to compute the vector
>
>F_i  =  - \int_Omega  \varphi_i(x)  Js(x)  dx
>
> for which you only need to be able to evaluate Js(x) and \varphi_i(x) at
> quadrature points x, but not a specific support points. This is no
> problem for the FE_NedelecSZ element. Can you explain why this does not
> work for you?
>
> 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.
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] Re: How to get the coordinates of a given degree of freedom on an edge?

2019-02-18 Thread Phạm Ngọc Kiên

Dear Prof. Wolfgang Bangerth,
I would like to thank you very much.
I am reading through the FE_NedelecSZ codes in Deal.II library.
I currently try to solve the electromagnetics problem, followed by the
research of Ross Kynch.
My goal is to compute the right hand side for the problem:
curl curl E + i*omega*mu*sigma*E = - Js
I only want to set the Js term here to be none zero, for example to be 1,
at a given degree of freedom (i.e. at source position), and  0 everywhere
else.
Thus, I need to write the codes to query if a degree of freedom is at
source position.
As you mentioned, now I think that I need to check if the source point is
on an edge when assembling the right hand side vector.
Could you please tell me help me about this?
Best regards,
Pham Ngoc Kien

Vào Th 3, 19 thg 2, 2019 vào lúc 15:45 Wolfgang Bangerth <
bange...@colostate.edu> đã viết:

> On 2/18/19 10:47 PM, Phạm Ngọc Kiên wrote:
> >
> > 2. I have tried some thing with the third method. And here-below is my
> code:
> >
> > Point p{0.5,0.5,0.5};//position in reference cell
> > Quadrature q(p);
> > FEValues fe_values_q(fe, q,update_quadrature_points);
> > fe_values_q.reinit(cell);
> >
> > //position in the real cell
> > std::vector> dof_pos = fe_values_q.get_quadrature_points();
> > for (int k =0; k  >  for (int j =0; j  >  std::cout< >  std::cout<<"\n";
> >  }
> > }
> >
> >  From this codes, I can get the real position of arbitrary points in my
> current reference cell.
> >
> > However, my purpose is to get the real position of a given degree of
> freedom.
> >
> > Could you please tell me how to get the position of a degree of freedom
> in the reference cell?
>
> The problem is that the NedelecSZ element does not have any such position.
> You
> are thinking of Lagrange elements (or the regular Nedelec elements) for
> which
> every shape function has a specific location where it is defined: It is
> one
> there, and zero at all other such locations. This is what we call "support
> points" in deal.II.
>
> But the NedelecSZ element is not defined this way (at least that's what I
> believe to be the case). Rather, shape functions are defined in such a way
> that certain *integrals* (rather than *point values*) are zero or one for
> all
> shape functions. As a consequence, there are no support points: Shape
> functions (and degrees of freedom) are not defined at individual points,
> but
> at whole edges.
>
> This means that your approach can not work. I don't know what you want to
> do,
> but let's start with that: What is your *goal*? Maybe then we can think
> about
> how to achieve it!
>
> 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.
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Re: [deal.II] Re: How to get the coordinates of a given degree of freedom on an edge?

2019-02-18 Thread Phạm Ngọc Kiên

I am sorry as my question is not clear enough.
What I mean is when we loop over all degree of freedom in a cell, for
example,

for (unsigned int i = 0; i < dofs_per_cell; ++i) { .. }

Do we have any method to get the position of the i_th degree of
freedom inside this loop?

I would like to thank you very much!

Best,

Pham Ngoc Kien


Vào Th 3, 19 thg 2, 2019 vào lúc 14:47 Phạm Ngọc Kiên <
ngockien.lea...@gmail.com> đã viết:

> Hi,
> I have tried with the methods described in the FAQ and figured out that:
> 1. As I used edge-based element (FE_NedelecSZ), I queried if my element
> had support points by  *fe.has_support_points()*. As a consequence, it
> returned zero.
> The *FiniteElement::get_unit_support_points() *also returned an empty
> vector.
> As there was no support point, DoFTools::map_dofs_to_support_points()
> could not work as you mentioned above.
>
> 2. I have tried some thing with the third method. And here-below is my
> code:
>
> Point p{0.5,0.5,0.5}; //position in reference cell
> Quadrature q(p);
> FEValues fe_values_q(fe, q, update_quadrature_points);
> fe_values_q.reinit(cell);
>
> //position in the real cell
> std::vector> dof_pos = fe_values_q.get_quadrature_points();
> for (int k = 0; k  for (int j = 0; j  std::cout< std::cout<<"\n";
> }
> }
>
> From this codes, I can get the real position of arbitrary points in my 
> current reference cell.
>
> However, my purpose is to get the real position of a given degree of freedom.
>
> Could you please tell me how to get the position of a degree of freedom in 
> the reference cell?
>
> Thank you very much for your help.
>
> Best,
>
> Pham Ngoc Kien
>
>
> Vào Th 3, 19 thg 2, 2019 vào lúc 11:47 Phạm Ngọc Kiên <
> ngockien.lea...@gmail.com> đã viết:
>
>> Yes, thank you very much.
>> I will dig deeper into this.
>> Best,
>> Pham Ngoc Kien
>>
>> Vào Th 3, 19 thg 2, 2019 vào lúc 11:23 Bruno Turcksin <
>> bruno.turck...@gmail.com> đã viết:
>>
>>> Pham Ngoc Kien,
>>>
>>> Le lun. 18 févr. 2019 à 19:33, Phạm Ngọc Kiên 
>>> a écrit :
>>>
>>>> With the DoFTools::map_dofs_to_support_points() , I saw in the library
>>>> manual that it can be used with no edge elements or the like.
>>>> As my finite element is edge-based one, I wonder if this function work?
>>>>
>>> No the function won't work but there are two other methods described in
>>> the FAQ that you could use. In particular, you can create a special
>>> quadrature at the points of interest.
>>>
>>> Best,
>>>
>>> Bruno
>>>
>>> --
>>> 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.
>>>
>>

-- 
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.

Re: [deal.II] Re: How to get the coordinates of a given degree of freedom on an edge?

2019-02-18 Thread Phạm Ngọc Kiên

Hi,
I have tried with the methods described in the FAQ and figured out that:
1. As I used edge-based element (FE_NedelecSZ), I queried if my element had
support points by  *fe.has_support_points()*. As a consequence, it returned
zero.
The *FiniteElement::get_unit_support_points() *also returned an empty
vector.
As there was no support point, DoFTools::map_dofs_to_support_points() could
not work as you mentioned above.

2. I have tried some thing with the third method. And here-below is my code:

Point p{0.5,0.5,0.5}; //position in reference cell
Quadrature q(p);
FEValues fe_values_q(fe, q, update_quadrature_points);
fe_values_q.reinit(cell);

//position in the real cell
std::vector> dof_pos = fe_values_q.get_quadrature_points();
for (int k = 0; k From this codes, I can get the real position of arbitrary points in my
current reference cell.

However, my purpose is to get the real position of a given degree of freedom.

Could you please tell me how to get the position of a degree of
freedom in the reference cell?

Thank you very much for your help.

Best,

Pham Ngoc Kien


Vào Th 3, 19 thg 2, 2019 vào lúc 11:47 Phạm Ngọc Kiên <
ngockien.lea...@gmail.com> đã viết:

> Yes, thank you very much.
> I will dig deeper into this.
> Best,
> Pham Ngoc Kien
>
> Vào Th 3, 19 thg 2, 2019 vào lúc 11:23 Bruno Turcksin <
> bruno.turck...@gmail.com> đã viết:
>
>> Pham Ngoc Kien,
>>
>> Le lun. 18 févr. 2019 à 19:33, Phạm Ngọc Kiên 
>> a écrit :
>>
>>> With the DoFTools::map_dofs_to_support_points() , I saw in the library
>>> manual that it can be used with no edge elements or the like.
>>> As my finite element is edge-based one, I wonder if this function work?
>>>
>> No the function won't work but there are two other methods described in
>> the FAQ that you could use. In particular, you can create a special
>> quadrature at the points of interest.
>>
>> Best,
>>
>> Bruno
>>
>> --
>> 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.
>>
>

-- 
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.

Re: [deal.II] Re: How to get the coordinates of a given degree of freedom on an edge?

2019-02-18 Thread Phạm Ngọc Kiên

Yes, thank you very much.
I will dig deeper into this.
Best,
Pham Ngoc Kien

Vào Th 3, 19 thg 2, 2019 vào lúc 11:23 Bruno Turcksin <
bruno.turck...@gmail.com> đã viết:

> Pham Ngoc Kien,
>
> Le lun. 18 févr. 2019 à 19:33, Phạm Ngọc Kiên 
> a écrit :
>
>> With the DoFTools::map_dofs_to_support_points() , I saw in the library
>> manual that it can be used with no edge elements or the like.
>> As my finite element is edge-based one, I wonder if this function work?
>>
> No the function won't work but there are two other methods described in
> the FAQ that you could use. In particular, you can create a special
> quadrature at the points of interest.
>
> Best,
>
> Bruno
>
> --
> 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.
>

-- 
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.

Re: [deal.II] Re: How to get the coordinates of a given degree of freedom on an edge?

2019-02-18 Thread Phạm Ngọc Kiên

Hi Bruno,
With the DoFTools::map_dofs_to_support_points() , I saw in the library
manual that it can be used with no edge elements or the like.
As my finite element is edge-based one, I wonder if this function work?
I have check the function fe.has_support_points() in my code and it
returned 0.
That is the reason why I put my question.
Thank you very much for your suggestion.

Best,
Pham Ngoc Kien

Vào Th 3, 19 thg 2, 2019 vào lúc 06:18 Bruno Turcksin <
bruno.turck...@gmail.com> đã viết:

> Hi,
>
> Take a look at this entry
> <https://github.com/dealii/dealii/wiki/Frequently-Asked-Questions#how-to-get-the-mapped-position-of-support-points-of-my-element>
> in the FAQ.
>
> Best,
>
> Bruno
>
> On Monday, February 18, 2019 at 12:20:48 AM UTC-5, Phạm Ngọc Kiên wrote:
>>
>> Dear colleagues,
>> I am working with FE_Nedelec_SZ element (edge-based finite element) in
>> Deal.II (version 9.0.0).
>> For my current problem, I need to get the global coordinates of a given
>> degree of freedom on and edge.
>> Do we have any code for this task?
>>
>> I have also thought of using material_id for an edge as an alternative
>> way to address my problem. However, I personally think that the material_id
>> can be used only with a cell.
>>
>> I would like to thank you very much in advance.
>>
>> Pham Ngoc Kien
>> Ph.D. student
>> Seoul National University
>>
> --
> 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.
>

-- 
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.

[deal.II] How to get the coordinates of a given degree of freedom on an edge?

2019-02-17 Thread Phạm Ngọc Kiên

Dear colleagues,
I am working with FE_Nedelec_SZ element (edge-based finite element) in 
Deal.II (version 9.0.0).
For my current problem, I need to get the global coordinates of a given 
degree of freedom on and edge. 
Do we have any code for this task?

I have also thought of using material_id for an edge as an alternative way 
to address my problem. However, I personally think that the material_id can 
be used only with a cell.

I would like to thank you very much in advance.

Pham Ngoc Kien
Ph.D. student
Seoul National University

-- 
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.

Re: [deal.II] output vector solution to vtk file

2019-01-25 Thread Phạm Ngọc Kiên

I really appreciate your answer. The library is really useful and i will
dig deeper into it.
I am using the latest version of dealii.
Thank you very much.
Best regards.
Pham Ngoc Kien

Vào 04:42 PM, T.6, 25 Th1, 2019 Jean-Paul Pelteret 
đã viết:

> Dear Pham,
>
> Welcome to deal.II!
>
> Your observation is related to this frequently asked question
> <https://github.com/dealii/dealii/wiki/Frequently-Asked-Questions#the-graphical-output-files-dont-make-sense-to-methey-seem-to-have-too-many-degrees-of-freedom>.
> Basically, deal.II outputs the solution on each cell individually because
> we cater for the situation where the solution is discontinuous between
> elements. This could arise when the solution field is discretised by
> discontinuous finite elements, when there are hanging nodes in the
> triangulation, or when a post-processed quantity is naturally discontinuous
> between elements. This patch-based output, although maybe a bit wasteful in
> the case of continuous fields due to the data duplication, is the most
> robust solution for all scenarios.
>
> For that particular file type and the default output options, the data
> points are located at the vertices of the triangulation. If you are using
> the developer version of deal.II it is now possible to output the solution
> to "higher-order cells” (e.g. quadratic or cubic interpolatory cells). In
> this case, you’d get extra data points on the edges and interior of each
> cell.
>
> I hope that this helped!
>
> Best,
> Jean-Paul
>
> On 25 Jan 2019, at 07:35, Phạm Ngọc Kiên 
> wrote:
>
> Dear Prof. Wolfgang Bangerth,
> I am a newbie with deal.II library, and I have a question in this topic.
> Thus, I put my question here.
> I can now output my vector values results using DataPostprocessor class
> into vtk file.
> However, when I opened the output file, I see some data points in the file
> had the same positions. I had checked with the examples (step-3) and
> received the same thing.
> Could you please tell me why it is?
> In my limited opinion, the function: attach_dof_handler(dof_handler) will
> give me the data positions which are the positions of all nodes in my mesh.
> I wonder if the same data points I mentioned above are the nodes on edge,
> cell, face?
> I would like to thank you very much in advance.
>
> Best regards,
> Pham Ngoc Kien
>
> --
> 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.
>
>
> --
> 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.
>

-- 
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.

Re: [deal.II] output vector solution to vtk file

2019-01-24 Thread Phạm Ngọc Kiên

Dear Prof. Wolfgang Bangerth,
I am a newbie with deal.II library, and I have a question in this topic. 
Thus, I put my question here.
I can now output my vector values results using DataPostprocessor class 
into vtk file.
However, when I opened the output file, I see some data points in the file 
had the same positions. I had checked with the examples (step-3) and 
received the same thing.
Could you please tell me why it is?
In my limited opinion, the function: attach_dof_handler(dof_handler) will 
give me the data positions which are the positions of all nodes in my mesh. 
I wonder if the same data points I mentioned above are the nodes on edge, 
cell, face?
I would like to thank you very much in advance.

Best regards,
Pham Ngoc Kien

-- 
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.

46 matches

Mail list logo