RE: [deal.II] Computing the solution gradient at the quadrature point on a face

2021-08-10 Thread Michael Li
Hi Jean-Pau, You’re always helpful. The deformation dependent loads are discussed in Wriggers, Wood and some other books and papers. I’m studying that. It seems the formulation is complicated. I believe that the auto-differentiation tools is a good approach to this nonlinear problem and I’ll dig into it later. For now I hope I can have a better understanding this topic and implement it in deal.II by hand then I will be comfortable with AD. Thank you!Michael  From: Jean-Paul PelteretSent: Sunday, August 8, 2021 12:16 PMTo: dealii@googlegroups.comSubject: Re: [deal.II] Computing the solution gradient at the quadrature point on a face Hi Michael, I’m glad that you’re on the right track now.  If you could give some advice or refer to some references on this topic, I would appreciate that very mush. Unfortunately, this is where pretty much I hit the limit of what I can say on this topic as it would conflict with the expectations of my employer. There’s a little by-line on one approach in this paperhttps://www.sciencedirect.com/science/article/abs/pii/S0997753814001508But as for other references, I’m hoping that someone else can give you some guidance there. (There’s at least one fairly common book in nonlinear mechanics that discusses the "linearisation of external virtual work”. Maybe a search in this direction might get you somewhere.) Another option would be to try the auto-differentiation tools to do all of the heavy lifting for you. With few restrictions, it would be compatible with any general traction load. At the very least and AD implementation it could help you validate any linearisation that you do by hand. You can take a look at steps 70 and 71 to get some idea as to what it might be able to do for you. Using the approach in step-71 is probably easiest to start off in your case. You could even limit it to just linearising the contribution that comes from this boundary / traction contribution. Best,Jean-Paul On 4. Aug 2021, at 23:22, Michael Li  wrote: Hi Jean-Paul, Yes, I set up the boundary check when taking care of the traction. It worked as expected using Physics::Transformations::nansons_formula(face_normal, F). To make a quick test, I did not linearize the load stiffness but just used the load at last Newton-Raphson step, so it is like a mixture of N-R and Picard iteration. The technique is not correct and there is some discrepancy compared with the reference. But the result is better than not considering the face area change at all. So in some sense, it tells me that the surface transform is correct. I’m working on how to linearize an arbitrary surface traction. I have found some references dealing with the normal(pressure) load but not arbitrary surface load (the force is not normal to the face in my case). If you could give some advice or refer to some references on this topic, I would appreciate that very mush. You’re totally right, it converges only if a very small load step is used. The residual grows up fast and the result is meaningless if the load increment is large. Thanks for your great help! Michael   From: Jean-Paul PelteretSent: Wednesday, August 4, 2021 12:19 AMTo: dealii@googlegroups.comSubject: Re: [deal.II] Computing the solution gradient at the quadrature point on a face Hi Michael, The one thing that stands out about the snippet of code that you shared is that there is no check to if (cell->face(face)->at_boundary() == true && )...Do you still have something like that? Otherwise, after a cursory glance, everything appears to be alright with this. It mimics what you’d do to get the deformation gradient at a QP in a cell.  To debug further, what happens when you apply a very small load? Maybe its not scaled appropriately for the stiffness of the material. Also, if you’ve now added a deformation dependent part to the load (i.e. the area transformation) then you might need to consistently linearise this in order to attain convergence for large loads. Best,Jean-PaulOn 2. Aug 2021, at 18:03, Michael Li  wrote: Hi Jean-Paul, Thank you for your hints.  I initialized a local variable solution_grads_u_face_total with fe_face_values_ref[u_fe].get_function_gradients to get the gradient of displacement at the quadrature point on the face (see codes below). But the gradient acquired is not sensible and the resultant deformation is totally wrong. void assemble_neumann_contribution_one_cell(…){    ….  const FEValuesExtractors::Vector _fe = data.solid->u_fe; std::vector2, dim,NumberType> > solution_grads_u_face_total(data.solid->qf_face.size()); scratch.fe_face_values_ref.reinit(cell, face); scratch.fe_face_values_ref[u_fe].get_function_gradients(scratch.solution_total,solution_grads_u_face_total);   for (unsigned int face = 0; face < GeometryInfo::faces_per_cell; ++face)   {   for (unsigned int f_q_point = 0; f_q_point < n_q_points_f; ++f_q_point)   {

[deal.II] Re: Boost version problem

2021-08-10 Thread Lucas Myers
Update: the problem is that, when I installed the newest version of dealii, 
the boost version that I had separately installed was 1.67.0. Instead of 
using the version bundled in dealii (version 1.70.0), dealii configured 
with the separately installed version. However, whenever it does the 
compilation check, it #includes the bundled version of boost. I'm not sure 
why it does that, because when I run a loose source file (e.g. 
`boost_test.cpp` in my home directory) with boost included, it uses the 
separately installed version. 

In any case, I think the fix will be to uninstall boost, and then 
re-install dealii so that it uses its bundled version. Just wanted to make 
a note about it because I think this is unintended behavior. 

- Lucas

On Tuesday, August 10, 2021 at 3:46:17 PM UTC-5 Lucas Myers wrote:

> Hi folks,
>
> I'm getting the error "The version number of boost that you are compiling 
> with does not match the version number of boost found during deal.II's 
> configuration step." when I try to compile the first tutorial program. It 
> is thrown at /usr/local/include/deal.II/base/config.h:508:17
>
> This is confusing to me, because when I go into 
> /usr/include/boost/version.hpp, it states that my BOOST_VERSION is 106700. 
> Further, in /usr/local/include/deal.II/base/config.h, it says that my 
> DEAL_II_BOOST_VERSION_MAJOR is 1, DEAL_II_BOOST_VERSION_MINOR is 67, and 
> DEAL_II_BOOST_VERSION_SUBMINOR is 0 which makes my 
> DEAL_II_BOOST_VERSION_GTE 106700, the same as what boost says it ought to 
> be.
>
> For context, I had installed dealii version 9.2.0 earlier, but had 
> problems when I tried to update it to include CUDA (I think there was some 
> API change with CUDA 11). I then installed dealii version 9.3.1 (without 
> doing anything to the other version). Now none of the example files in 
> version 9.3.1 will compile.
>
> Any help is appreciated, thanks so much
>
> Kind regards,
> Lucas
>

-- 
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/ff8355c1-75c9-4d0c-bf90-014357e1ac37n%40googlegroups.com.


[deal.II] Boost version problem

2021-08-10 Thread Lucas Myers
Hi folks,

I'm getting the error "The version number of boost that you are compiling 
with does not match the version number of boost found during deal.II's 
configuration step." when I try to compile the first tutorial program. It 
is thrown at /usr/local/include/deal.II/base/config.h:508:17

This is confusing to me, because when I go into 
/usr/include/boost/version.hpp, it states that my BOOST_VERSION is 106700. 
Further, in /usr/local/include/deal.II/base/config.h, it says that my 
DEAL_II_BOOST_VERSION_MAJOR is 1, DEAL_II_BOOST_VERSION_MINOR is 67, and 
DEAL_II_BOOST_VERSION_SUBMINOR is 0 which makes my 
DEAL_II_BOOST_VERSION_GTE 106700, the same as what boost says it ought to 
be.

For context, I had installed dealii version 9.2.0 earlier, but had problems 
when I tried to update it to include CUDA (I think there was some API 
change with CUDA 11). I then installed dealii version 9.3.1 (without doing 
anything to the other version). Now none of the example files in version 
9.3.1 will compile.

Any help is appreciated, thanks so much

Kind regards,
Lucas

-- 
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/04760d59-d6ef-464a-8f5d-d5b63cd7e31fn%40googlegroups.com.


[deal.II] finite element with shape functions having delta_ij property at quadrature points

2021-08-10 Thread Simon

Dear all,

I am solving a problem in 2d using FE_Q(2) elements and a gauss quadrature 
rule with (fe.degree +1) quadrature points in each co-ordinate direction, 
that is, I have in total nine quadrature points. My question pertains to 
the following:
At each cell, I need to approximate a field whose sampling (support) points 
are the quadrature points belonging to reduced integration, i.e, there are 
four quadrature points in my case and the four (shape) functions 
approximating my field should be designed as follows:
N_j (xi_k) = delta_{jk} ,
where xi_k are the coordinates of the quadrature points. So I need four 
(shape) functions, each of which is one at one of the four quadrature 
points and zero at the three others.  
That said, my ansatz is given by (the coefficents a(xi) are of course known)
f(xi) = a(xi_1) N_1(xi) + a(xi_2) N_2(xi) + a(xi_3) N_3(xi) + a(xi_4) N(xi)
and I need to evaluate the function f(xi) at the *nine* quadrature points.

What is the best way to do set up the FEValues object?
I have seen that there is a constructor for the FE_Q element which takes a 
Quadratute<1> object. I guess this would help me to define the (shape) 
functions pertaining to the field f(xi), but I think I can not evaluate 
this field at the nine quadrature points, because (i) their local 
coordinates are of course different in the new FEValues object and (ii) 
second, I would have to insert negative local coordinates for a set of 
them. 
Maybe I do not even need a FEValues object for my purpuse. As I said, I 
only need to do the approximation f(xi) and evaluate it at the nine 
quadrature points.

Any input is greatly appreciated!

Best
Simon



-- 
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/852af19c-dd9c-4c94-809f-e7b8b6df49d9n%40googlegroups.com.


Re: [deal.II] point_value, Real parts, step-29

2021-08-10 Thread Jean-Paul Pelteret
Another thing: You probably also need to initialise with the right number of 
components. So something like
Vector vecSol (fe.n_components());


> On 10. Aug 2021, at 19:40, Jean-Paul Pelteret  wrote:
> 
> Hi Hermes,
> 
> You don’t say what errors you’re seeing, but my guess is that it now doesn’t 
> compile. This variant of the function (the one that Daniel linked to) returns 
> void, so you should call it before outputting the result:
> 
> Vector vecSol;
> VectorTools::point_value(dof_handler, solution,Point<2>(0.2, 0.2),vecSol)
> std::cout << "Solution at (0.2,0.2): "<< vecSol << std::endl;
> 
> This will hopefully get you what you want.
> 
> Best,
> Jean-Paul
> 
> 
>> On 10. Aug 2021, at 16:09, Hermes Sampedro > > wrote:
>> 
>> Thank you for your answer. I am not sure if I fully understand your 
>> suggestion. Do you mean something like that:
>> 
>>  Vector vecSol;
>>  std::cout << "Solution at (0.2,0.2): "<< 
>> VectorTools::point_value(dof_handler, solution,Point<2>(0.2, 0.2),vecSol)<< 
>> std::endl;
>>
>> 
>> I still get some errors. Is there not any way to get for example the real 
>> part of solution easily and use it directly on the point_value as in step-3?
>> 
>> Thank you 
>> H
>> 
>> El mar, 10 ago 2021 a las 15:49, Daniel Arndt (> >) escribió:
>> Hermes,
>> 
>> Use another overload. The one returning the solution as a parameter should 
>> work: 
>> https://www.dealii.org/current/doxygen/deal.II/namespaceVectorTools.html#acd358e9b110ccbf4a7f76796d206b9c7
>>  
>> 
>> 
>> Best,
>> Daniel
>> 
>> Am Di., 10. Aug. 2021 um 09:41 Uhr schrieb Hermes Sampedro 
>> mailto:hermesampe...@gmail.com>>:
>> Dear all, 
>> 
>> It is explained in Step-3 how to evaluate the solution in a point. I am 
>> trying to do the same for Step-29, to evaluate the real and imaginary parts 
>> separately in a single point:
>> std::cout << "Solution at (0.2,0.2): "<< 
>> VectorTools::point_value(dof_handler, solution,Point<2>(0.2, 0.2))<< 
>> std::endl;
>> 
>> I do not have any problems compiling however, an error occurs when running:
>> 
>> The violated condition was:  dof.get_fe(0).n_components() == 1
>> 
>> What is the proper way to call the real and imaginary parts of the solution 
>> at a particular point here?
>> 
>> 
>> 
>> Thank you very much!
>> 
>> H.
>> 
>> 
>> -- 
>> 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/54a2d44e-194a-43c0-aa7f-9fddd5bd0dean%40googlegroups.com
>>  
>> .
>> 
>> -- 
>> 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 a topic in the 
>> Google Groups "deal.II User Group" group.
>> To unsubscribe from this topic, visit 
>> https://groups.google.com/d/topic/dealii/Ru1_uMbix30/unsubscribe 
>> .
>> To unsubscribe from this group and all its topics, send an email to 
>> dealii+unsubscr...@googlegroups.com 
>> .
>> To view this discussion on the web visit 
>> https://groups.google.com/d/msgid/dealii/CAOYDWb%2BR0bRTVxB-F22Hwk_ZexN4%3DVwCbZ7ZKqs3JjEaHFhncA%40mail.gmail.com
>>  
>> .
>> 
>> -- 
>> 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 
>> 

Re: [deal.II] point_value, Real parts, step-29

2021-08-10 Thread Jean-Paul Pelteret
Hi Hermes,

You don’t say what errors you’re seeing, but my guess is that it now doesn’t 
compile. This variant of the function (the one that Daniel linked to) returns 
void, so you should call it before outputting the result:

Vector vecSol;
VectorTools::point_value(dof_handler, solution,Point<2>(0.2, 0.2),vecSol)
std::cout << "Solution at (0.2,0.2): "<< vecSol << std::endl;

This will hopefully get you what you want.

Best,
Jean-Paul


> On 10. Aug 2021, at 16:09, Hermes Sampedro  wrote:
> 
> Thank you for your answer. I am not sure if I fully understand your 
> suggestion. Do you mean something like that:
> 
>  Vector vecSol;
>  std::cout << "Solution at (0.2,0.2): "<< 
> VectorTools::point_value(dof_handler, solution,Point<2>(0.2, 0.2),vecSol)<< 
> std::endl;
>
> 
> I still get some errors. Is there not any way to get for example the real 
> part of solution easily and use it directly on the point_value as in step-3?
> 
> Thank you 
> H
> 
> El mar, 10 ago 2021 a las 15:49, Daniel Arndt ( >) escribió:
> Hermes,
> 
> Use another overload. The one returning the solution as a parameter should 
> work: 
> https://www.dealii.org/current/doxygen/deal.II/namespaceVectorTools.html#acd358e9b110ccbf4a7f76796d206b9c7
>  
> 
> 
> Best,
> Daniel
> 
> Am Di., 10. Aug. 2021 um 09:41 Uhr schrieb Hermes Sampedro 
> mailto:hermesampe...@gmail.com>>:
> Dear all, 
> 
> It is explained in Step-3 how to evaluate the solution in a point. I am 
> trying to do the same for Step-29, to evaluate the real and imaginary parts 
> separately in a single point:
> std::cout << "Solution at (0.2,0.2): "<< 
> VectorTools::point_value(dof_handler, solution,Point<2>(0.2, 0.2))<< 
> std::endl;
> 
> I do not have any problems compiling however, an error occurs when running:
> 
> The violated condition was:  dof.get_fe(0).n_components() == 1
> 
> What is the proper way to call the real and imaginary parts of the solution 
> at a particular point here?
> 
> 
> 
> Thank you very much!
> 
> H.
> 
> 
> -- 
> 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/54a2d44e-194a-43c0-aa7f-9fddd5bd0dean%40googlegroups.com
>  
> .
> 
> -- 
> 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 a topic in the Google 
> Groups "deal.II User Group" group.
> To unsubscribe from this topic, visit 
> https://groups.google.com/d/topic/dealii/Ru1_uMbix30/unsubscribe 
> .
> To unsubscribe from this group and all its topics, send an email to 
> dealii+unsubscr...@googlegroups.com 
> .
> To view this discussion on the web visit 
> https://groups.google.com/d/msgid/dealii/CAOYDWb%2BR0bRTVxB-F22Hwk_ZexN4%3DVwCbZ7ZKqs3JjEaHFhncA%40mail.gmail.com
>  
> .
> 
> -- 
> 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/CAB%3DnHhZrnXUeYt5FrKa66oTwoV8WJ8b%3DFqLFJ9pfP%3DUHRTyV%2BQ%40mail.gmail.com
>  
> .

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

[deal.II] deal.II Newsletter #177

2021-08-10 Thread 'Rene Gassmoeller' via deal.II User Group
Hello everyone!

This is deal.II newsletter #177.
It automatically reports recently merged features and discussions about the 
deal.II finite element library.


## Below you find a list of recently proposed or merged features:

#12639: Replace AssertDimension (proposed by gfcas) 
https://github.com/dealii/dealii/pull/12639

#12638: Integrate difference (proposed by gfcas) 
https://github.com/dealii/dealii/pull/12638

#12637: Typos (proposed by gfcas; merged) 
https://github.com/dealii/dealii/pull/12637

#12636: step-40 disable timer summary output (proposed by gfcas) 
https://github.com/dealii/dealii/pull/12636

#12635: Override MGTransferGlobalCoarsening::prolongate_and_add() (proposed by 
peterrum; merged) https://github.com/dealii/dealii/pull/12635

#12634: Fix MGTransferMatrixFree Assert [WIP] (proposed by tjhei) 
https://github.com/dealii/dealii/pull/12634

#12632: Step-47: unify code (proposed by peterrum; merged) 
https://github.com/dealii/dealii/pull/12632

#12631: Small modifications of FEEvaluation documentation (proposed by 
peterrum; merged) https://github.com/dealii/dealii/pull/12631

#12630: Add face and cell queries to FEFaceValuesBase and FEInterfaceValues 
(proposed by jppelteret; merged) https://github.com/dealii/dealii/pull/12630

#12629: Fix signal residual_step (proposed by peterrum; merged) 
https://github.com/dealii/dealii/pull/12629

#12628: Add dof_indices() to FEInterfaceValues (proposed by jppelteret; merged) 
https://github.com/dealii/dealii/pull/12628

#12627: Update test output for CUDA (proposed by masterleinad; merged) 
https://github.com/dealii/dealii/pull/12627

#12625: Integrate hessians on cells (proposed by bergbauer; merged) 
https://github.com/dealii/dealii/pull/12625

#12624: Avoid a warning about copying invalid data. (proposed by bangerth) 
https://github.com/dealii/dealii/pull/12624

#12623: Added diagnostic functions for `AffineConstraints`. (proposed by 
marcfehling; merged) https://github.com/dealii/dealii/pull/12623

#12622: Add another example to the ScalarFunctionFromFunctionObject 
documentation. (proposed by bangerth) 
https://github.com/dealii/dealii/pull/12622

#12618: Extend documentation of TriangulationDescription (proposed by peterrum; 
merged) https://github.com/dealii/dealii/pull/12618

#12616: Switch from unsigned int to ConstraintTypes (proposed by peterrum; 
merged) https://github.com/dealii/dealii/pull/12616

#12615: Extend test (proposed by peterrum; merged) 
https://github.com/dealii/dealii/pull/12615

#12614: Standardized output for many `mpi/hp_unify_dof_indices_X` tests. 
(proposed by marcfehling; merged) https://github.com/dealii/dealii/pull/12614

#12593: Evaluate and integrate hessians (proposed by bergbauer; merged) 
https://github.com/dealii/dealii/pull/12593

#12411: FEInterfaceValues: Make class interfaces consistent with FEValues and 
its views (proposed by jppelteret; merged) 
https://github.com/dealii/dealii/pull/12411


## And this is a list of recently opened or closed discussions:

#12633: MGTransferMatrixFree Assert failure (opened) 
https://github.com/dealii/dealii/issues/12633

#12626: Missing identity constraints in 3D. (opened) 
https://github.com/dealii/dealii/issues/12626

#12621: symmetric_tensor.templates.h (opened) 
https://github.com/dealii/dealii/issues/12621

#12117: deal.II Release 9.3.0 (closed) 
https://github.com/dealii/dealii/issues/12117


A list of all major changes since the last release can be found at 
https://www.dealii.org/developer/doxygen/deal.II/recent_changes.html.


Thanks for being part of the community!


Let us know about questions, problems, bugs or just share your experience by 
writing to dealii@googlegroups.com, or by opening issues or pull requests at 
https://www.github.com/dealii/dealii.
Additional information can be found at https://www.dealii.org/.

-- 
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/6112a2a0.1c69fb81.30c66.e3d1SMTPIN_ADDED_MISSING%40gmr-mx.google.com.


Re: [deal.II] point_value, Real parts, step-29

2021-08-10 Thread Hermes Sampedro
Thank you for your answer. I am not sure if I fully understand your
suggestion. Do you mean something like that:

 Vector<*double*> vecSol;

 std::cout << "Solution at (0.2,0.2): "<<
VectorTools::point_value(dof_handler, solution,Point<2>(0.2, 0.2),vecSol)<<
std::endl;




I still get some errors. Is there not any way to get for example the real
part of *solution* easily and use it directly on the point_value as in
step-3?


Thank you

H

El mar, 10 ago 2021 a las 15:49, Daniel Arndt ()
escribió:

> Hermes,
>
> Use another overload. The one returning the solution as a parameter should
> work:
> https://www.dealii.org/current/doxygen/deal.II/namespaceVectorTools.html#acd358e9b110ccbf4a7f76796d206b9c7
>
> Best,
> Daniel
>
> Am Di., 10. Aug. 2021 um 09:41 Uhr schrieb Hermes Sampedro <
> hermesampe...@gmail.com>:
>
>> Dear all,
>>
>> It is explained in Step-3 how to evaluate the solution in a point. I am
>> trying to do the same for Step-29, to evaluate the real and imaginary parts
>> separately in a single point:
>>
>> *std::cout << "Solution at (0.2,0.2): "<<
>> VectorTools::point_value(dof_handler, solution,Point<2>(0.2, 0.2))<<
>> std::endl;*
>>
>> I do not have any problems compiling however, an error occurs when
>> running:
>>
>> *The violated condition was:  dof.get_fe(0).n_components() == 1*
>>
>> What is the proper way to call the real and imaginary parts of the
>> solution at a particular point here?
>>
>>
>> Thank you very much!
>>
>> H.
>>
>> --
>> 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/54a2d44e-194a-43c0-aa7f-9fddd5bd0dean%40googlegroups.com
>> 
>> .
>>
> --
> 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 a topic in the
> Google Groups "deal.II User Group" group.
> To unsubscribe from this topic, visit
> https://groups.google.com/d/topic/dealii/Ru1_uMbix30/unsubscribe.
> To unsubscribe from this group and all its topics, send an email to
> dealii+unsubscr...@googlegroups.com.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/dealii/CAOYDWb%2BR0bRTVxB-F22Hwk_ZexN4%3DVwCbZ7ZKqs3JjEaHFhncA%40mail.gmail.com
> 
> .
>

-- 
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/CAB%3DnHhZrnXUeYt5FrKa66oTwoV8WJ8b%3DFqLFJ9pfP%3DUHRTyV%2BQ%40mail.gmail.com.


Re: [deal.II] point_value, Real parts, step-29

2021-08-10 Thread Daniel Arndt
Hermes,

Use another overload. The one returning the solution as a parameter should
work:
https://www.dealii.org/current/doxygen/deal.II/namespaceVectorTools.html#acd358e9b110ccbf4a7f76796d206b9c7

Best,
Daniel

Am Di., 10. Aug. 2021 um 09:41 Uhr schrieb Hermes Sampedro <
hermesampe...@gmail.com>:

> Dear all,
>
> It is explained in Step-3 how to evaluate the solution in a point. I am
> trying to do the same for Step-29, to evaluate the real and imaginary parts
> separately in a single point:
>
> *std::cout << "Solution at (0.2,0.2): "<<
> VectorTools::point_value(dof_handler, solution,Point<2>(0.2, 0.2))<<
> std::endl;*
>
> I do not have any problems compiling however, an error occurs when running:
>
> *The violated condition was:  dof.get_fe(0).n_components() == 1*
>
> What is the proper way to call the real and imaginary parts of the
> solution at a particular point here?
>
>
> Thank you very much!
>
> H.
>
> --
> 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/54a2d44e-194a-43c0-aa7f-9fddd5bd0dean%40googlegroups.com
> 
> .
>

-- 
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/CAOYDWb%2BR0bRTVxB-F22Hwk_ZexN4%3DVwCbZ7ZKqs3JjEaHFhncA%40mail.gmail.com.


[deal.II] point_value, Real parts, step-29

2021-08-10 Thread Hermes Sampedro
Dear all, 

It is explained in Step-3 how to evaluate the solution in a point. I am 
trying to do the same for Step-29, to evaluate the real and imaginary parts 
separately in a single point:

*std::cout << "Solution at (0.2,0.2): "<< 
VectorTools::point_value(dof_handler, solution,Point<2>(0.2, 0.2))<< 
std::endl;*

I do not have any problems compiling however, an error occurs when running:

*The violated condition was:  dof.get_fe(0).n_components() == 1*

What is the proper way to call the real and imaginary parts of the solution 
at a particular point here?


Thank you very much!

H.

-- 
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/54a2d44e-194a-43c0-aa7f-9fddd5bd0dean%40googlegroups.com.


Re: [deal.II] Inaccurate convergence rate using ParsedConvergenceTable class

2021-08-10 Thread luca.heltai
If you want to use symbolic calculations, you could also leverage SymEngine, 
and use Functions::SymbolicFunction

https://www.dealii.org/current/doxygen/deal.II/classFunctions_1_1SymbolicFunction.html

In your code, you could create such object by making a 

std::unique_ptr> fun;

parse its expression from file, and then 

fun = std::make_unique>(expression);

In this way the generated function would compute the gradient symbolically.

L.

> On 9 Aug 2021, at 16:39, Marco Feder  wrote:
> 
> You're right, I totally missed the default parameter h in the constructor.
> 
> Thanks,
> Marco
> Il giorno domenica 8 agosto 2021 alle 22:26:01 UTC+2 Timo Heister ha scritto:
> FunctionParser uses a finite difference approximation for the gradient. I 
> think this is explained in the class documentation.
> This explains the results you see...
> 
> On Sun, Aug 8, 2021, 15:15 Marco Feder  wrote:
> Ciao Luca!
> 
> > Can you check if you get the same results with this order?
> Yes, the rates are the same. 
> 
> > Are you calling `error_from_exact` with mapping as a first argument?
> Yes
> 
> Actually, I've fixed the issue after reading this:
> > If the Solution class implements the Gradient, then 
> > ParsedConvergenceTable should use that. 
> 
> Indeed, I realised I wrote:
> error_table.error_from_exact(mapping, dof_handler, solution, exact_solution);
> 
> where exact_solution is a FunctionParser which will be initialised with 
> the expression of the analytical solution, so it doesn't give any info about 
> the Gradient. Instead, if I replace exact_solution in the last argument with 
> the "classical" Solution() (so that there's also an implementation of 
> Gradient) I have my expected rates. Apparently it wasn't using the Gradient 
> function. 
> 
> However, looking at the source code and in the documentation of 
> integrate_difference, I don't see how the H^1 norm was computed without 
> having an implementation of the Gradient. I would expect to see an exception 
> asking for its implementation. Was it using some FD-like approximation or am 
> I missing something?
> 
> 
> Thanks,
> Marco
> 
> 
> Il giorno domenica 8 agosto 2021 alle 18:52:03 UTC+2 luca@gmail.com ha 
> scritto:
> The other option is that you are not calling the `error_from_exact` function 
> that takes a mapping argument, but the one without it, and your boundary 
> cells are introducing some error due to a worse mapping…. Are you calling 
> `error_from_exact` with mapping as a first argument?
> 
> Luca
> 
>> Il giorno 8 ago 2021, alle ore 6:49 PM, Luca Heltai  ha 
>> scritto:
>> 
>> 
>> Ciao Marco,
>> 
>>> I was expanding step-12 using a manufactured solution to check the order p 
>>> in H^1 (and p+1 in L^2) norm on a uniformly refined mesh for the DG upwind 
>>> method. I've used the ParsedConvergenceTable with a parameter file as 
>>> described in the docs. As Rate_key I'm using the DoFs, while as Rate_mode I 
>>> have reduction_rate_log2. 
>>> 
>>> With p=1 and p=2 everything is fine. But if I set the finite element degree 
>>> to 3, then the H^1 convergence rate decreases, as you can see in the 
>>> attached image.
>>> 
>> 
>> Is it possible that you are using different quadrature rules in the two 
>> cases? Your image shows a deterioration of the error on the order of 1e-8 
>> for H1, and 1e-12 for L2, which is very close to machine precision.
>> 
>> Internally the parsed convergence table does the exact same thing you wrote 
>> explicitly. (If you check the source code, you’ll see that it calls 
>> integrate_difference and compute_global_error for each error type you 
>> specify in the parameter file).
>> 
>>> This, however, doesn't happen if I use a classical ConvergenceTable. 
>>> Namely, I first compute the local error in each cell, and then the global 
>>> error in the classical way:
>>> VectorTools::integrate_difference(mapping,dof_handler,solution, 
>>> Solution(),H1_error_per_cell, 
>>> QGauss(fe->tensor_degree()+1),VectorTools::H1_norm);
>>> 
>>> const double H1_error = VectorTools::compute_global_error(triangulation, 
>>> H1_error_per_cell,  VectorTools::H1_norm); //assuming I provided also the 
>>> gradient method for the Solution class
>>> 
>>> 
>>> 
>>> Does anyone have any idea why this is happening? My guess was that while 
>>> computing the H^1 semi-norm the ParsedConvergenceTable class does some 
>>> approximation to compute the gradient from the exact solution expression 
>>> and hence that could be the source of the issue. Conversely, in the 
>>> "ConvergenceTable" way I do define explicitly the gradient of the exact 
>>> solution in the Solution class.
>> 
>> If the Solution class implements the Gradient, then 
>> ParsedConvergenceTable should use that. 
>> 
>> You are calling “error_from_exact” of that class, right?
>> 
>> This is where it is called:
>> 
>> https://www.dealii.org/current/doxygen/deal.II/parsed__convergence__table_8h_source.html
>> 
>> Line 621. As you see, the quadrature rule used is a Gauss formula of