Re: [deal.II] try catch of "ExcDereferenceInvalidCell"

[Simon Wiesheier]

" So I am not sure how you would ever use an `Assert(...)` in optimized
release mode (where you want to run your program to get any
performance). "

I was indeed considering my program in debug mode.
Based on what you said, there is no way for me to catch
the failed .reinit call in optimized released mode
and using if-statements is the proper way for me
to handle my conditions.
Is that correct?

" Have you considered looking into our FEPointEvaluation [1] and related
"nonmatching" facilities?"

Yes, I already use FEPointEvaluation to evaluate the interpoland
at arbitrary points.
What other "nonmatching" facilities do you think may be helpful?
As said, I have to evaluate the interpoland (first and second derivatives)
at arbitrary points in real space --
as often as there are quadrature points on my geometry in each assembly.

Best,
Simon



Am Fr., 4. Aug. 2023 um 00:02 Uhr schrieb Matthias Maier :

> Have you considered looking into our FEPointEvaluation [1] and related
> "nonmatching" facilities?
>
> I have the feeling that we already provide the right data structures for
> your interpolation/evaluation needs - in particular with an efficient
> and fast implementation. An example usage can be found in step-19.
>
> Best,
> Matthias
>
> [1]
> https://www.dealii.org/current/doxygen/deal.II/classFEPointEvaluation.html
>
> --
> 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/87zg37bx9e.fsf%4043-1.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/CAM50jEuBtewGX962LmiPtQyLPbz%2BScvDJ_7akdtFQFdmmVcoOQ%40mail.gmail.com.

Re: [deal.II] try catch of "ExcDereferenceInvalidCell"

[Matthias Maier]

Have you considered looking into our FEPointEvaluation [1] and related
"nonmatching" facilities?

I have the feeling that we already provide the right data structures for
your interpolation/evaluation needs - in particular with an efficient
and fast implementation. An example usage can be found in step-19.

Best,
Matthias

[1]  https://www.dealii.org/current/doxygen/deal.II/classFEPointEvaluation.html

-- 
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/87zg37bx9e.fsf%4043-1.org.

Re: [deal.II] try catch of "ExcDereferenceInvalidCell"

[Matthias Maier]

I am a bit confused.


In optimized release mode the `Assert` are reduced to a nnop and will
not abort (or throw an exception). For user-visible/induced runtime
errors we have `AssertThrow` which will throw an exception (and not
abort a program). In summary:

 - we use `Assert(...)` in the library to check for logical error
   conditions in *debug* mode. These are optimized out in optimized
   release mode.

   The reasoning behind this is that these asserts are "static", logical
   mistakes in the library / user code.

 - we use `AssertThrow(...)` to check for error conditions in *debug*
   and *release* mode. These checks are always run and used to validate
   user-controlled input.

   (and catchable to do something in case such an error is encountered).

So I am not sure how you would ever use an `Assert(...)` in optimized
release mode (where you want to run your program to get any
performance).


Regarding the additional if statement: If this is performance critical I
suggest to not use the c++ exception mechanism in this way. Stack
unwinding is very, very slow.


On the other hand, if you have a very unlikely error path and handled in
an if-statement (c++20 style):

if (/* error condition */) [[unlikely]] {
  // ... handle error
}

then we have often made the observation that this leads to very acceptable
performance - even in hot execution paths.


Best,
Matthias




On Thu, Aug  3, 2023, at 12:25 CDT, Simon Wiesheier  
wrote:

> Let me please clarify my last concerns in this regard:
>
> The operations
> GridTools::find_active_cell_around_point and
> feValues.reinit
> are executed at quadrature point level and performance
> is really important in our application.
>
> The reason why I wanted to use an
> try/catch block is to surrogate the use of an if-statement.
> Consider these two variants:
>
> Variant 1:
>
> disable_abort_on_exception();
> try
> {
>feValues.reinit(cell,...);
> }
> catch(dealii::Exception Base exception & )
> {
>Assert(checkSecondRunTimeCondition, ...)
>// do something
> }
>
> Variant 2:
>
> if(cell->state()==-1)
> {
> Assert(checkSecondRunTimeCondition, ...)
>// do something
> }
> else
> {
>  feValues.reinit(cell,...);
> }
>
> As you can see, there are two conditions I have to check at every
> quadrature point.
> Clearly, Variant 2 implements an if-requests and is probably not the way to
> go.
>
>
> 1. So with regard to performance, would you prefer Variant 1
> (as the try/catch introduces less overhead compared to the if-statement)?
> 2. If so, is the call of "disable_abort_on_exception()" good
> coding practice or would you refrain from doing so?
> In the main function of our program, we have a try/catch with
> catch(dealii::ExceptionBase exception &).
> Given that, in my opinion  "disable_abort_on_exception()"
> should cause any undue behavior.
>
> Best,
> Simon
>
>
> Am Do., 3. Aug. 2023 um 14:21 Uhr schrieb Wolfgang Bangerth <
> bange...@colostate.edu>:
>
>> On 8/2/23 15:58, Simon Wiesheier wrote:
>> >
>> > Does that make sense?
>> > Or do you see a better solution, like checking
>> > cell->status() right after GridTools::find_active_cell_around_point?
>>
>> This. If the function you call returns an error code in the form of an end
>> iterator, just test for that rather than doing something with this
>> iterator
>> and wait to see what happens. It is always worth checking errors as early
>> as
>> possible.
>>
>> 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/53307567-6cd6-3f21-2c63-a5eedf5367f5%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/874jlfdc1x.fsf%4043-1.org.

[deal.II] How to compute derivatives of the composite function with the reference cell mapping

[Anton]

Hello,

I would like to use the NonMatching::QuadratureGenerator class 
directly, as I need to generate multiple (as in many!, unfortunately) 
nonmatching quadratures for the same cell.  Therefore I would like to avoid 
the route of interpolating the level set function onto the whole mesh etc, 
as described in CutFEM tutorial.

Therefore I am thinking about calling the "generate" function, which only 
needs a level set function and a bounding box as inputs.  According to the 
documentation the level set function should provide the gradient and the 
Hessian.  I would like to do this in the reference cell coordinates, so 
that I can simply use this quadrature later via FEView as one normally does.
 
 * The bounding box I presume is the hypercube [0,1]^{dim} for the meshes I 
am interested in. (Can one request them for the reference cell and not the 
triangulation cell?)

 * The level set function should also be reasonably straightforward.  I 
know the function in the physical coordinates (say, f(x)) and all its 
derivatives.  I can also get the map from reference to physical coordinates 
via FEView, let us call this map x=M(y).  The issue is how to extract the 
derivatives of this map so that the chain rule can be applied?  I am mostly 
working with Q1 maps, so of course I could figure out the derivatives "by 
hand".  I am just wondering if there is a more intelligent (code-wise) way 
of doing this?

-- 
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/4ac38563-e467-4f1d-8cb4-ee35379515a6n%40googlegroups.com.

Re: [deal.II] Neumann Boundary condition in MeshWorker for Automatic Differentiation

[Abbas Ballout]

Jost, 

For looping on the boundaries, you'll have to write an extra lambda amd 
pass that to the MeshWorker. 
You can look at step 12 or 74 for this.
Also following this cause I want to know what are the MeshWorker 
alternatives since I am relying heavily on it in my code

Best, 
Abbas
 
On Thursday, August 3, 2023 at 3:11:48 PM UTC+2 jost@googlemail.com 
wrote:

> Thanks for the quick response!
>
> Would you suggest me to read more tutorials on how the MeshWorker works to 
> fix this issue then or try to set up my program without a MeshWorker? I 
> thought I had set it up correctly and tried already for a while to figure 
> out what I did wrong. Also I Initially only started using the MeshWorker 
> because it seemed to be recommended in Step-72 and I thought following the 
> recommendations is probably the best and easiest way. Now I am not so sure 
> anymore about it.
>
> Best,
>
> Jost
>
> Wolfgang Bangerth schrieb am Donnerstag, 3. August 2023 um 14:49:44 UTC+2:
>
>> On 8/3/23 04:54, 'Jost Arndt' via deal.II User Group wrote: 
>> > 
>> > for (unsigned int f = 0; f < GeometryInfo::faces_per_cell; ++f) 
>> > { 
>> > if 
>> > (cell->face(f)->at_boundary())//(face->at_boundary()) 
>> > { 
>> > std::cout << "reached line 
>> 1032" << 
>> > std::endl; 
>> > const auto _face_values = 
>> > scratch_data.reinit(cell, f); 
>> > std::cout << "reached line 
>> 1034" << 
>> > std::endl; 
>> >  /* here integration using 
>> > fe_face_value should happen*/ 
>> > 
>> > but the last cout never gets executed, but the code crashes. I get a 
>> long but 
>> > cryptic (to me?) stacktrace but also the following error message: 
>> > 
>> > 
>> > An error occurred in line <3044> of file <./source/fe/fe_values.cc> in 
>> function 
>> > dealii::FEValuesBase::FEValuesBase(unsigned int, 
>> unsigned int, 
>> > dealii::UpdateFlags, const dealii::Mapping&, const 
>> > dealii::FiniteElement&) [with int dim = 2; int spacedim 
>> >  = 2] 
>> > The violated condition was: 
>> > n_q_points > 0 
>> > Additional information: 
>> > There is nothing useful you can do with an FEValues object when 
>> using 
>> > a quadrature formula with zero quadrature points! 
>> > 
>> > From here I am quite unsure what I did wrong or how else to create the 
>> > FEFaceValues Object? I would be really thankful about any hint! 
>>
>> The error likely means that MeshWorker hasn't initialized the 
>> FEFaceValues 
>> object. Nobody knows any more how MeshWorker really works, and so we try 
>> to 
>> discourage use of that framework. But I do know that MeshWorker allows 
>> you to 
>> not only loop over the cells, but also the faces of a triangulation. That 
>> is 
>> the recommended way to deal with face integrals (rather than dealing with 
>> these face integrals as part of cell integration -- because you want to 
>> make 
>> sure that you visit each face only once, rather than twice in your 
>> approach) 
>> and if you tell MeshWorker that you also want to use a face worker, then 
>> it 
>> should (?) initialize the FEFaceValues object as well. 
>>
>> Best 
>> W. 
>>
>> -- 
>>  
>> Wolfgang Bangerth email: bang...@colostate.edu 
>> www: http://www.math.colostate.edu/~bangerth/ 
>>
>>
>>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/90721349-725a-49e7-976c-7d989be0812bn%40googlegroups.com.

Re: [deal.II] try catch of "ExcDereferenceInvalidCell"

[Simon Wiesheier]

Let me please clarify my last concerns in this regard:

The operations
GridTools::find_active_cell_around_point and
feValues.reinit
are executed at quadrature point level and performance
is really important in our application.

The reason why I wanted to use an
try/catch block is to surrogate the use of an if-statement.
Consider these two variants:

Variant 1:

disable_abort_on_exception();
try
{
   feValues.reinit(cell,...);
}
catch(dealii::Exception Base exception & )
{
   Assert(checkSecondRunTimeCondition, ...)
   // do something
}

Variant 2:

if(cell->state()==-1)
{
Assert(checkSecondRunTimeCondition, ...)
   // do something
}
else
{
 feValues.reinit(cell,...);
}

As you can see, there are two conditions I have to check at every
quadrature point.
Clearly, Variant 2 implements an if-requests and is probably not the way to
go.


1. So with regard to performance, would you prefer Variant 1
(as the try/catch introduces less overhead compared to the if-statement)?
2. If so, is the call of "disable_abort_on_exception()" good
coding practice or would you refrain from doing so?
In the main function of our program, we have a try/catch with
catch(dealii::ExceptionBase exception &).
Given that, in my opinion  "disable_abort_on_exception()"
should cause any undue behavior.

Best,
Simon


Am Do., 3. Aug. 2023 um 14:21 Uhr schrieb Wolfgang Bangerth <
bange...@colostate.edu>:

> On 8/2/23 15:58, Simon Wiesheier wrote:
> >
> > Does that make sense?
> > Or do you see a better solution, like checking
> > cell->status() right after GridTools::find_active_cell_around_point?
>
> This. If the function you call returns an error code in the form of an end
> iterator, just test for that rather than doing something with this
> iterator
> and wait to see what happens. It is always worth checking errors as early
> as
> possible.
>
> 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/53307567-6cd6-3f21-2c63-a5eedf5367f5%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/CAM50jEsCfErRV8-Su0gw4MyJS13W%3DwfbY5mCZY2_LXYMcm8efw%40mail.gmail.com.

Re: [deal.II] Neumann Boundary condition in MeshWorker for Automatic Differentiation

['Jost Arndt' via deal.II User Group]

Thanks for the quick response!

Would you suggest me to read more tutorials on how the MeshWorker works to 
fix this issue then or try to set up my program without a MeshWorker? I 
thought I had set it up correctly and tried already for a while to figure 
out what I did wrong. Also I Initially only started using the MeshWorker 
because it seemed to be recommended in Step-72 and I thought following the 
recommendations is probably the best and easiest way. Now I am not so sure 
anymore about it.

Best,

Jost

Wolfgang Bangerth schrieb am Donnerstag, 3. August 2023 um 14:49:44 UTC+2:

> On 8/3/23 04:54, 'Jost Arndt' via deal.II User Group wrote:
> > 
> > for (unsigned int f = 0; f < GeometryInfo::faces_per_cell; ++f)
> > {
> > if 
> > (cell->face(f)->at_boundary())//(face->at_boundary())
> > {
> > std::cout << "reached line 1032" 
> << 
> > std::endl;
> > const auto _face_values = 
> > scratch_data.reinit(cell, f);
> > std::cout << "reached line 1034" 
> << 
> > std::endl;
> >  /* here integration using 
> > fe_face_value should happen*/
> > 
> > but the last cout never gets executed, but the code crashes. I get a 
> long but 
> > cryptic (to me?) stacktrace but also the following error message:
> > 
> > 
> > An error occurred in line <3044> of file <./source/fe/fe_values.cc> in 
> function
> > dealii::FEValuesBase::FEValuesBase(unsigned int, unsigned 
> int, 
> > dealii::UpdateFlags, const dealii::Mapping&, const 
> > dealii::FiniteElement&) [with int dim = 2; int spacedim
> >  = 2]
> > The violated condition was:
> > n_q_points > 0
> > Additional information:
> > There is nothing useful you can do with an FEValues object when using
> > a quadrature formula with zero quadrature points!
> > 
> > From here I am quite unsure what I did wrong or how else to create the 
> > FEFaceValues Object? I would be really thankful about any hint!
>
> The error likely means that MeshWorker hasn't initialized the FEFaceValues 
> object. Nobody knows any more how MeshWorker really works, and so we try 
> to 
> discourage use of that framework. But I do know that MeshWorker allows you 
> to 
> not only loop over the cells, but also the faces of a triangulation. That 
> is 
> the recommended way to deal with face integrals (rather than dealing with 
> these face integrals as part of cell integration -- because you want to 
> make 
> sure that you visit each face only once, rather than twice in your 
> approach) 
> and if you tell MeshWorker that you also want to use a face worker, then 
> it 
> should (?) initialize the FEFaceValues object as well.
>
> Best
> W.
>
> -- 
> 
> Wolfgang Bangerth email: bang...@colostate.edu
> www: http://www.math.colostate.edu/~bangerth/
>
>
>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/cb043712-8b46-4b60-9abe-09dc25f6383en%40googlegroups.com.

Re: [deal.II] Neumann Boundary condition in MeshWorker for Automatic Differentiation

[Wolfgang Bangerth]

On 8/3/23 04:54, 'Jost Arndt' via deal.II User Group wrote:


for (unsigned int f = 0; f < GeometryInfo::faces_per_cell; ++f)
                         {
                                 if 
(cell->face(f)->at_boundary())//(face->at_boundary())

                                 {
                                         std::cout << "reached line 1032" << 
std::endl;
                                         const auto _face_values = 
scratch_data.reinit(cell, f);
                                         std::cout << "reached line 1034" << 
std::endl;
  /* here integration using 
fe_face_value should happen*/


but the last cout never gets executed, but the code crashes. I get a long but 
cryptic (to me?) stacktrace but also the following error message:



An error occurred in line <3044> of file <./source/fe/fe_values.cc> in function
dealii::FEValuesBase::FEValuesBase(unsigned int, unsigned int, 
dealii::UpdateFlags, const dealii::Mapping&, const 
dealii::FiniteElement&) [with int dim = 2; int spacedim

  = 2]
The violated condition was:
     n_q_points > 0
Additional information:
     There is nothing useful you can do with an FEValues object when using
     a quadrature formula with zero quadrature points!

 From here I am quite unsure what I did wrong or how else to create the 
FEFaceValues Object? I would be really thankful about any hint!


The error likely means that MeshWorker hasn't initialized the FEFaceValues 
object. Nobody knows any more how MeshWorker really works, and so we try to 
discourage use of that framework. But I do know that MeshWorker allows you to 
not only loop over the cells, but also the faces of a triangulation. That is 
the recommended way to deal with face integrals (rather than dealing with 
these face integrals as part of cell integration -- because you want to make 
sure that you visit each face only once, rather than twice in your approach) 
and if you tell MeshWorker that you also want to use a face worker, then it 
should (?) initialize the FEFaceValues object as well.


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/4b0ffb76-3000-788f-6f90-2b98f3c94b0f%40colostate.edu.

Re: [deal.II] Q: Question about extracting part of a vector

[Wolfgang Bangerth]

On 8/3/23 05:01, Najwa Alshehri wrote:
The first scenario involves directly using the FeFieldFunction when my 
computed solution consists of a single component. In this case, I have 
successfully calculated the L2 and H1 norms of the error using the 
`integrate_difference` and `compute_global_error` functions. The calculations 
proceeded smoothly, and the execution speed was reasonable.


However, the second scenario involves using the FeFieldFunction when the 
computed solution consists of two components. My objective is to calculate the 
error for the first component only. To accomplish this, I employed the new 
VectorFunctionFromScalarFunctionObject class locally. Unfortunately, the 
computation of the L2 norm in this case was considerably slow and the H1 norm 
was very slow.


Given the current situation, I have two questions that I hope you can assist 
me with, as you have done in the past:


1. Regarding my specific project, is there an alternative method, other than 
the FeFieldFunction, to save the deserialized solution as a function and 
utilize it as the exact solution? It is crucial that the solution be saved in 
the first component of a vector with two components.


2. In a more general sense, how can I test this function to determine if the 
sluggishness is solely attributable to the FeFieldFunction? Additionally, I 
would like to verify whether the function is functioning correctly or not, to 
be able to send it in a pull request to the dealii git repository.


Najwa:
I don't know why the run time should be substantially slower for a system of 
two components than for a single scalar problem. That might be worth 
investigating. If you wanted to look into this, it would be useful to start 
with as small a program that shows the issue.


In general, FEFieldFunction is slow because it doesn't know where the function 
should be evaluated. If you give it a point, it first has to find the cell the 
point is in, and then compute the reference coordinates of that point, etc. 
You probably do this many times if you want to evaluate the solution on the 
quadrature points of a whole triangulation. This will never be fast


There are perhaps better ways to deal with the situation. For example, instead 
of doing the integration on the coarse mesh, you could take the coarse 
solution, interpolate it onto the fine mesh, and then do the integration over 
the fine mesh -- because there you can then just compute the values of both 
the coarse and fine solution at the quadrature points of the fine mesh using 
FEValues. This is fast, and the interpolation to the fine mesh should also be 
reasonably fast.


There are also faster but more specialized alternatives to FEFieldFunction. 
FEPointEvaluation is one option, though I don't know that class well enough.


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/e79e6d07-c8ea-170a-abd9-8e0cedf6ca9d%40colostate.edu.

Re: [deal.II] try catch of "ExcDereferenceInvalidCell"

[Wolfgang Bangerth]

On 8/2/23 15:58, Simon Wiesheier wrote:


Does that make sense?
Or do you see a better solution, like checking
cell->status() right after GridTools::find_active_cell_around_point?


This. If the function you call returns an error code in the form of an end 
iterator, just test for that rather than doing something with this iterator 
and wait to see what happens. It is always worth checking errors as early as 
possible.


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/53307567-6cd6-3f21-2c63-a5eedf5367f5%40colostate.edu.

Re: [deal.II] Q: Question about extracting part of a vector

[Najwa Alshehri]

Dear team,

 

I want to express my gratitude for the support you have provided thus far. 
I would like to provide you with an update on the current situation and 
seek your assistance regarding a couple of issues.

 

Firstly, I have successfully implemented the missing function 
`VectorFunctionFromScalarFunctionObject::gradient()` in my project. 
Although it appears to be functioning correctly, the execution speed is 
extremely slow. In fact, it took until the next working day to obtain the 
results for only two cycles of refinement. Based on the response from 
@Wolfgang, it seems that the issue lies with the FeFieldFunction.

 

Allow me to provide some context. I have a reference solution serialized on 
a fine mesh. In my current code, my objective is to deserialize this 
solution and utilize it as the exact solution. To achieve this, I have 
employed the FeFieldFunction. Now, I have two scenarios in which I need to 
utilize the FeFieldFunction.

 

The first scenario involves directly using the FeFieldFunction when my 
computed solution consists of a single component. In this case, I have 
successfully calculated the L2 and H1 norms of the error using the 
`integrate_difference` and `compute_global_error` functions. The 
calculations proceeded smoothly, and the execution speed was reasonable.

 

However, the second scenario involves using the FeFieldFunction when the 
computed solution consists of two components. My objective is to calculate 
the error for the first component only. To accomplish this, I employed the 
new VectorFunctionFromScalarFunctionObject class locally. Unfortunately, 
the computation of the L2 norm in this case was considerably slow and the 
H1 norm was very slow.

 

Given the current situation, I have two questions that I hope you can 
assist me with, as you have done in the past:

 

1. Regarding my specific project, is there an alternative method, other 
than the FeFieldFunction, to save the deserialized solution as a function 
and utilize it as the exact solution? It is crucial that the solution be 
saved in the first component of a vector with two components.

 

2. In a more general sense, how can I test this function to determine if 
the sluggishness is solely attributable to the FeFieldFunction? 
Additionally, I would like to verify whether the function is functioning 
correctly or not, to be able to send it in a pull request to the dealii git 
repository.

 

Your guidance and expertise in resolving these issues would be greatly 
appreciated.

 

Thank you once again for your continued support.

 

Best regards,

Najwa

On Tuesday, August 1, 2023 at 12:50:27 PM UTC+3 abbas.b...@gmail.com wrote:

> Najwa, 
> The files I edited are in the pull request here: 
> https://github.com/dealii/dealii/pull/15805.
> Let me know if I can help in a better way or of there is something I can 
> do. 
>
> Abbas. 
>
> On Tuesday, August 1, 2023 at 5:48:41 AM UTC+2 najwaa...@gmail.com wrote:
>
>> Hello Abbas,
>>
>> Thank you for your response. I have actually written something, but I 
>> need to test it. It would be great to have a look at your work as well, so 
>> we can compare.
>>
>> Best regards,
>> Najwa
>>
>> On Mon, 31 Jul 2023 at 10:48 PM Abbas Ballout  
>> wrote:
>>
>>> Najwaa, 
>>>
>>> I submitted a pull request  
>>> recently about something similar 
>>> . (I guess)
>>> Maybe that would help. 
>>>
>>> Abbas 
>>>  
>>>
>>> On Monday, July 31, 2023 at 7:37:59 PM UTC+2 najwaa...@gmail.com wrote:
>>>
 Dear Wolfgang,

 Thank you for your answer. I will try to write a patch to the deal.II 
 sources that implement the missing function. This might require some time. 
 I will come back here if needed.


 Thank you all for your quick answers,
 Najwa

 On Wednesday, July 26, 2023 at 11:41:15 PM UTC+3 Wolfgang Bangerth 
 wrote:

> On 7/26/23 12:17, Najwa Alshehri wrote: 
> > 
> > Regarding the function "VectorFucntionFromScalarFunctionObject," I 
> have 
> > observed that it works for computing the L2 norm of the error. 
> However, when I 
> > attempted to compute the H1_seminorm of the error, I noticed that 
> the gradient 
> > was not implemented for this function. Please let me know if my 
> understanding 
> > is incorrect, as I referred to the  for 
> this. 
>
> Yes, correct. Someone needs to implement 
> VectorFucntionFromScalarFunctionObject::gradient() 
> and/or 
> VectorFucntionFromScalarFunctionObject::gradient_list() 
> in exactly the same way as the value() and value_list() functions are 
> implemented. 
>
>
> > Could you kindly advise me on the simplest way to accomplish this? 
>
> The easiest way is to create a class derived from 
> VectorFucntionFromScalarFunctionObject in your own project that 
> implements the 
> missing 

[deal.II] Neumann Boundary condition in MeshWorker for Automatic Differentiation

['Jost Arndt' via deal.II User Group]

Hi there,

so far I had a great learning curve, but now i seem to be stuck.
I try to solve my own little PDE, which is a system of multiple nonlinear 
time dependent equations and read through plenty of the Tutorials. 
Using Dirichlet-Boundary conditions everything worked out finally, now I 
wanted to switch to Neuman boundary conditions and I get some issue I do 
not understand.

I followed mostly step-72 to understand how to use automatic 
differentiation. For Neumann boundary conditions I have now to integrate 
over the boundary of my domain, therefore I jumped right into the 
cell_worker lambda function (step-72) right after the for-loop over all 
quadrature points. 

There I tried to stick to the example given in here: 
https://dealii.org/developer/doxygen/deal.II/classMeshWorker_1_1ScratchData.html
 
and created a for loop:

for (unsigned int f = 0; f < GeometryInfo::faces_per_cell; ++f) 
{
if 
(cell->face(f)->at_boundary())//(face->at_boundary())
{
std::cout << "reached line 1032" << 
std::endl;
const auto _face_values = 
scratch_data.reinit(cell, f);
std::cout << "reached line 1034" << 
std::endl;
 /* here integration using 
fe_face_value should happen*/

but the last cout never gets executed, but the code crashes. I get a long 
but cryptic (to me?) stacktrace but also the following error message:


An error occurred in line <3044> of file <./source/fe/fe_values.cc> in 
function   

dealii::FEValuesBase::FEValuesBase(unsigned int, unsigned 
int, dealii::UpdateFlags, const dealii::Mapping&, const 
dealii::FiniteElement&) [with int dim = 2; int spacedim
 = 2]   

  
The violated condition was: 

  
n_q_points > 0 

   
Additional information: 

  
There is nothing useful you can do with an FEValues object when using   

  
a quadrature formula with zero quadrature points!

>From here I am quite unsure what I did wrong or how else to create the 
FEFaceValues Object? I would be really thankful about any hint! 

-- 
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/03776b07-f004-4099-9428-36feb17d706fn%40googlegroups.com.