Re: [deal.II] Facilities for computing the "goodness" of a mesh (i.e. something like a mesh-quality indicator)

2020-01-15 Thread Wolfgang Bangerth 
On 1/15/20 10:42 AM, Krishnakumar Gopalakrishnan wrote:
> 
> Thank you. This is certainly looking useful.  Do we have any example usage of 
> this function? For example, how do we re-compute the mesh based on the return 
> value from this aspect_ratio function or something?

Other than the functions Jean-Paul referred to, and the 
GridTools::laplace_transform() function, there is nothing in deal.II for this 
purpose.

But I took the liberty to write a few statements in the discussion in step-6 
on this topic:
   https://github.com/dealii/dealii/pull/9344

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/0dfcbdfe-0d80-5492-0a39-72ef57c8db7b%40colostate.edu.

Re: [deal.II] Using solution on from one FE problem as a boundary condition for another

2020-01-15 Thread Jean-Paul Pelteret 
Hey Ernesto,

I’ve used the approach outlined by Wolfgang and Daniel in my own work (multi 
physics problems, single DoFHandler with BC’s being transferred between a 
primary and secondary problem solved on subsets of the domain). Overall, it 
works really well so I’d back up their suggestion and recommend that you adopt 
it. Its mostly really just a question of finding common support points on the 
two surfaces, and seeing which DoF indices link up and transferring data from 
one solution vector to the other. I was looking though my code to see if I 
could share a snippet with you, but its so close in structure to the 
interpolate_boundary_values() (and still retains all of the original comments 
in the deal.II source :-) that its not really worth it. 

Cheers,
J-P

> On 15 Jan 2020, at 07:48, Daniel Arndt  wrote:
> 
> Ernesto,
> 
> Imposing strong boundary conditions should not be overly difficult in this 
> case. 
> Just as Wolgang said you should have a look at interpolate_boundary_values() 
> and replace the Function object evaluations on each cell by the function 
> values given through
> a FEValues object initialized for your first solution. Assuming that you are 
> using interpolatory ansatz functions, this interpolation has to happen at the 
> support points of your target finite element.
> In particular, these support points should be the quadrature points for the 
> FEValues object.
> 
> Best,
> Daniel
> 
> Am Di., 14. Jan. 2020 um 01:50 Uhr schrieb Ernesto Ismail 
> mailto:ernesto.ism...@gmail.com>>:
> Dear Prof. Bangerth,
> 
> Thank you for your reply. 
> 
> I would ideally have liked to impose the boundary conditions strongly but it 
> looks like doing so will not be straight-forward. I will do some 
> investigations into applying the boundary conditions weakly and attempt that 
> first.
> 
> One approach I had tried was to use Lobatto quadrature and use a quadrature 
> point history data structure to transfer the solution at the boundary. The 
> problem with this approach was that the differing element types lead to some 
> strange behaviours along the boundary.
> 
> Thanks again,
> Ernesto
> 
> On Thursday, 9 January 2020 23:43:27 UTC+2, Wolfgang Bangerth wrote:
> 
> Ernesto, 
> 
> > I am currently working with a staggered coupled system where I am solving 
> > several finite element problems in sequence. Up until now my boundary 
> > conditions have been fairly straightforward and I've been able to get 
> > meaningful results from my work. 
> > 
> > The finite element schemes I use for each step are fairly different, with 
> > some 
> > being discontinuous and of varying orders. 
> > 
> > I am now in a position where I need to use the solution of one of my steps 
> > as 
> > the boundary condition for the next (I actually need to combine two of them 
> > with a simple arithmetic function, but baby steps for now). I presume I 
> > need 
> > to use some sort of projection to achieve this, but I'm not exactly sure 
> > how 
> > to go about it. I thought to try and use the project_boundary_values() 
> > function to do what I need, but I am unsure how I would structure 
> > the boundary_functions argument. 
> 
> Do you need to enforce these boundary conditions strongly? If you can enforce 
> them weakly, then the solution of one problem only enters the weak form of 
> another problem (in that case through a boundary integral) and that is really 
> not very different from the way we do this in other problems where we have 
> multiple DoFHandlers -- e.g., in step-31. (Though my usual advice still holds 
> true that I think everything becomes simpler if you pack everything into one 
> DoFHandler). 
> 
> It gets more complicated if you want to impose strong boundary conditions. 
> You've already found the FEFieldFunction function, but it's complicated and 
> slow. A simpler way would be if you looked at the implementation of the 
> interpolate_boundary_values() function and identified where it is asking the 
> Function object for its values. That's the place where you'd need to somehow 
> combine your existing solutions into something else and turn that into a 
> constraint. 
> 
> 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

Re: [deal.II] Facilities for computing the "goodness" of a mesh (i.e. something like a mesh-quality indicator)

2020-01-15 Thread Jean-Paul Pelteret 
Dear Krishna,

There is also this long-standing pull request
https://github.com/dealii/dealii/pull/8196 

that hooks into the Verdict library to do mesh quality assessment. One day, 
hopefully in the not too distant future, I’ll get back to integrating it into 
the deal.II library. But you could always extract the salient portions of that 
code for your purposes.

Best,
Jean-Paul

> On 15 Jan 2020, at 18:31, Martin Kronbichler  
> wrote:
> 
> Dear Krishna,
> 
> We have recently added one such function, namely 
> GridTools::compute_aspect_ratio_of_cells(),
> 
> https://dealii.org/developer/doxygen/deal.II/namespaceGridTools.html#a2a9d8801d952ea8afa9ec1573229a329
>  
> 
> Best,
> Martin
> On 15.01.20 18:26, Krishnakumar Gopalakrishnan wrote:
>> In the results/extension section of Step-6 as well in Step-49, there is a 
>> great deal of discussion about obtaining good meshes are for a problem, with 
>> the caveat that this is a highly specialized topic. Step-49 is full of 
>> warnings about things that can go wrong when applying certain 
>> transformations, and it seems to implant doubts in my mind that my meshes 
>> are never going to be correct :)
>> 
>> 
>> Given this, in the spirit of the Kelly Indicator for solution errors in a 
>> cell, is there a similar single-metric for the quality/goodness of each cell 
>> in a mesh? Maybe there is something based on the smallest internal angle (to 
>> prevent slivers, for Jacobian of the mapping to be invertible, to prevent 
>> inverted cells etc).
>> 
>> 
>> After all, error estimation is a complicated topic. nevertheless the Kelly 
>> Indicator seems to do a reasonable job for quite a few equations. If 
>> currently not implemented in Deal.II, is there any journal article that 
>> discusses research in this direction?
>> 
>> 
>> 
>> Krishna
>> 
>> -- 
>> 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/3602a53a-a98a-4a5c-a512-a6d179c8dcec%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/36574d22-6fa9-e153-637a-2641e798756d%40gmail.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/2138DAC3-68CF-46D7-9E2F-12DA93702150%40gmail.com.

Re: [deal.II] discontinous contour over elements

2020-01-15 Thread David Eaton 
Bruno,

Thanks for your explanation. Indeed, the gradient of velocity, vorticity, will 
be discontinuous over C0 elements. I will give a try and do a L2 projection. 
Surprisingly I did not experience it in my code before. I used standard 
Lagrangian shape functions and continuous Galerkin method. I suppose dealii 
used the same interpolation functions. I have to investigate why my simple code 
does not have this discontinuous vorticity field.

Thanks
D.


From: dealii@googlegroups.com  on behalf of Bruno 
Blais 
Sent: Thursday, January 16, 2020, 1:07 AM
To: deal.II User Group
Subject: Re: [deal.II] discontinous contour over elements

An easy way to carry the projection that Wolfgang suggested is to use an L2 
projection.
The L2 projection matrix is only a mass matrix and your right hand side is 
constructed by the integral of multiplication of the variable you want to 
project with the test function. Generally, this matrix is very very easy to 
invert. This will yield you the C0 representation of your discontinuous field 
(vorticity) such that the error between your C0 projection and your original 
field is minimized at the position of the gauss points.

This is a procedure we use to set the initial conditions in Lethe when the 
initial condition is a complex function (for instance a Taylor-Green vortex).


On Wednesday, 15 January 2020 11:22:52 UTC-5, David Eaton wrote:
I understand the C0 element is piecewise linear across elements. However, I did 
not experience the same issue in my own C++ code while I use C0 element with 
the Petrov Galerkin stabilization terms. Actually, I am very confused at this 
point. How could I get rid of it while using C0 element?

Thanks
D.


From: dea...@googlegroups.com  on behalf of Wolfgang 
Bangerth 
Sent: Thursday, January 16, 2020, 12:05 AM
To: dea...@googlegroups.com
Subject: Re: [deal.II] discontinous contour over elements

On 1/14/20 10:04 PM, David Eaton wrote:
>
> Thank you for your suggestions. I am going to take a look at Lethe and compare
> with my implementation. In stabilized formulation, I used quadrilateral
> element, instead of P2 P1 Taylor-Hood element. The used element is only C0
> element. I also did not expect such a discontinuity between elements. Although
> I use P2 P1 Taylor-Hood element without stabilization terms, the discontinuity
> is still there. Probably I made mistakes somewhere  in setup. I also suspect
> that my solution is not converged. After taking a small relative tolerance
> 1e-8, the discontinuity still appears.

David, you did not understand what we were saying: If you use C0 elements
(think, piecewise linear) and you take derivatives to compute the vorticity,
then you automatically get a discontinuous function. That has nothing to do
with stabilization, solver tolerances, etc. It's just a consequence of the
fact that C0 elements and their shape functions have kinks and consequently
their derivatives are discontinuous.

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 dea...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/f03c64fc-3736-340a-1cdc-f24749fb413f%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/bbe79d0e-d423-4df8-8f83-e4e934fedd7c%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/SL2PR01MB2553926E27A0DC8D3638BCECC3370%40SL2PR01MB2553.apcprd01.prod.exchangelabs.com.

Re: [deal.II] discontinous contour over elements

2020-01-15 Thread David Eaton 
Professor Bangerth,

Thanks for the hints. I almost got what you mean. I will try to perform a L2 
projection.

Regards
David

From: dealii@googlegroups.com  on behalf of Wolfgang 
Bangerth 
Sent: Thursday, January 16, 2020, 1:03 AM
To: dealii@googlegroups.com
Subject: Re: [deal.II] discontinous contour over elements

On 1/15/20 9:22 AM, David Eaton wrote:
> I understand the C0 element is piecewise linear across elements. However, I
> did not experience the same issue in my own C++ code while I use C0 element
> with the Petrov Galerkin stabilization terms. Actually, I am very confused at
> this point. How could I get rid of it while using C0 element?

Not easily. But you can project the discontinuous function (namely, the
gradient or vorticity of your velocity) onto a continuous finite element
field. Then you again have a continuous function -- albeit a different one, of
course, from the one you actually wanted.

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/4c39585b-55dd-c495-17b3-950b1e0f79d5%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/SL2PR01MB255349EE978C91452641B41CC3370%40SL2PR01MB2553.apcprd01.prod.exchangelabs.com.

Re: [deal.II] Facilities for computing the "goodness" of a mesh (i.e. something like a mesh-quality indicator)

2020-01-15 Thread Krishnakumar Gopalakrishnan 
Hi Martin,

Thank you. This is certainly looking useful.  Do we have any example usage 
of this function? For example, how do we re-compute the mesh based on the 
return value from this aspect_ratio function or something?

Krishna

On Wednesday, January 15, 2020 at 5:31:34 PM UTC, Martin Kronbichler wrote:
>
> Dear Krishna,
>
> We have recently added one such function, namely 
> GridTools::compute_aspect_ratio_of_cells(),
>
>
> https://dealii.org/developer/doxygen/deal.II/namespaceGridTools.html#a2a9d8801d952ea8afa9ec1573229a329
>
> Best,
> Martin
> On 15.01.20 18:26, Krishnakumar Gopalakrishnan wrote:
>
> In the results/extension section of Step-6 as well in Step-49, there is a 
> great deal of discussion about obtaining good meshes are for a problem, 
> with the caveat that this is a highly specialized topic. Step-49 is full of 
> warnings about things that can go wrong when applying certain 
> transformations, and it seems to implant doubts in my mind that my meshes 
> are never going to be correct :)
>
>
> Given this, in the spirit of the Kelly Indicator for solution errors in a 
> cell, is there a similar single-metric for the quality/goodness of each 
> cell in a mesh? Maybe there is something based on the smallest internal 
> angle (to prevent slivers, for Jacobian of the mapping to be invertible, to 
> prevent inverted cells etc).
>
>
> After all, error estimation is a complicated topic. nevertheless the Kelly 
> Indicator seems to do a reasonable job for quite a few equations. If 
> currently not implemented in Deal.II, is there any journal article that 
> discusses research in this direction?
>
>
>
> Krishna
> -- 
> 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 dea...@googlegroups.com .
> To view this discussion on the web visit 
> https://groups.google.com/d/msgid/dealii/3602a53a-a98a-4a5c-a512-a6d179c8dcec%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/0d406456-3ef9-4624-adb5-0343602466f4%40googlegroups.com.

Re: [deal.II] Facilities for computing the "goodness" of a mesh (i.e. something like a mesh-quality indicator)

2020-01-15 Thread Martin Kronbichler 

Dear Krishna,

We have recently added one such function, namely 
GridTools::compute_aspect_ratio_of_cells(),


https://dealii.org/developer/doxygen/deal.II/namespaceGridTools.html#a2a9d8801d952ea8afa9ec1573229a329

Best,
Martin

On 15.01.20 18:26, Krishnakumar Gopalakrishnan wrote:


In the results/extension section of Step-6 as well in Step-49, there 
is a great deal of discussion about obtaining good meshes are for a 
problem, with the caveat that this is a highly specialized topic. 
Step-49 is full of warnings about things that can go wrong when 
applying certain transformations, and it seems to implant doubts in my 
mind that my meshes are never going to be correct :)



Given this, in the spirit of the Kelly Indicator for solution errors 
in a cell, is there a similar single-metric for the quality/goodness 
of each cell in a mesh? Maybe there is something based on the smallest 
internal angle (to prevent slivers, for Jacobian of the mapping to be 
invertible, to prevent inverted cells etc).



After all, error estimation is a complicated topic. nevertheless the 
Kelly Indicator seems to do a reasonable job for quite a few 
equations. If currently not implemented in Deal.II, is there any 
journal article that discusses research in this direction?




Krishna

--
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/3602a53a-a98a-4a5c-a512-a6d179c8dcec%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/36574d22-6fa9-e153-637a-2641e798756d%40gmail.com.

[deal.II] Facilities for computing the "goodness" of a mesh (i.e. something like a mesh-quality indicator)

2020-01-15 Thread Krishnakumar Gopalakrishnan 


In the results/extension section of Step-6 as well in Step-49, there is a 
great deal of discussion about obtaining good meshes are for a problem, 
with the caveat that this is a highly specialized topic. Step-49 is full of 
warnings about things that can go wrong when applying certain 
transformations, and it seems to implant doubts in my mind that my meshes 
are never going to be correct :)


Given this, in the spirit of the Kelly Indicator for solution errors in a 
cell, is there a similar single-metric for the quality/goodness of each 
cell in a mesh? Maybe there is something based on the smallest internal 
angle (to prevent slivers, for Jacobian of the mapping to be invertible, to 
prevent inverted cells etc).


After all, error estimation is a complicated topic. nevertheless the Kelly 
Indicator seems to do a reasonable job for quite a few equations. If 
currently not implemented in Deal.II, is there any journal article that 
discusses research in this direction?



Krishna

-- 
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/3602a53a-a98a-4a5c-a512-a6d179c8dcec%40googlegroups.com.

Re: [deal.II] discontinous contour over elements

2020-01-15 Thread Bruno Blais 
An easy way to carry the projection that Wolfgang suggested is to use an L2 
projection.
The L2 projection matrix is only a mass matrix and your right hand side is 
constructed by the integral of multiplication of the variable you want to 
project with the test function. Generally, this matrix is very very easy to 
invert. This will yield you the C0 representation of your discontinuous 
field (vorticity) such that the error between your C0 projection and your 
original field is minimized at the position of the gauss points.

This is a procedure we use to set the initial conditions in Lethe when the 
initial condition is a complex function (for instance a Taylor-Green 
vortex).


On Wednesday, 15 January 2020 11:22:52 UTC-5, David Eaton wrote:
>
> I understand the C0 element is piecewise linear across elements. However, 
> I did not experience the same issue in my own C++ code while I use C0 
> element with the Petrov Galerkin stabilization terms. Actually, I am very 
> confused at this point. How could I get rid of it while using C0 element?
>
> Thanks
> D.
>
> --
> *From:* dea...@googlegroups.com   > on behalf of Wolfgang Bangerth  >
> *Sent:* Thursday, January 16, 2020, 12:05 AM
> *To:* dea...@googlegroups.com 
> *Subject:* Re: [deal.II] discontinous contour over elements
>
> On 1/14/20 10:04 PM, David Eaton wrote:
> > 
> > Thank you for your suggestions. I am going to take a look at Lethe and 
> compare 
> > with my implementation. In stabilized formulation, I used quadrilateral 
> > element, instead of P2 P1 Taylor-Hood element. The used element is only 
> C0 
> > element. I also did not expect such a discontinuity between elements. 
> Although 
> > I use P2 P1 Taylor-Hood element without stabilization terms, the 
> discontinuity 
> > is still there. Probably I made mistakes somewhere  in setup. I also 
> suspect 
> > that my solution is not converged. After taking a small relative 
> tolerance 
> > 1e-8, the discontinuity still appears.
>
> David, you did not understand what we were saying: If you use C0 elements 
> (think, piecewise linear) and you take derivatives to compute the 
> vorticity, 
> then you automatically get a discontinuous function. That has nothing to 
> do 
> with stabilization, solver tolerances, etc. It's just a consequence of the 
> fact that C0 elements and their shape functions have kinks and 
> consequently 
> their derivatives are discontinuous.
>
> 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 dea...@googlegroups.com .
> To view this discussion on the web visit 
> https://groups.google.com/d/msgid/dealii/f03c64fc-3736-340a-1cdc-f24749fb413f%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/bbe79d0e-d423-4df8-8f83-e4e934fedd7c%40googlegroups.com.

Re: [deal.II] discontinous contour over elements

2020-01-15 Thread Wolfgang Bangerth 
On 1/15/20 9:22 AM, David Eaton wrote:
> I understand the C0 element is piecewise linear across elements. However, I 
> did not experience the same issue in my own C++ code while I use C0 element 
> with the Petrov Galerkin stabilization terms. Actually, I am very confused at 
> this point. How could I get rid of it while using C0 element?

Not easily. But you can project the discontinuous function (namely, the 
gradient or vorticity of your velocity) onto a continuous finite element 
field. Then you again have a continuous function -- albeit a different one, of 
course, from the one you actually wanted.

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/4c39585b-55dd-c495-17b3-950b1e0f79d5%40colostate.edu.

Re: [deal.II] discontinous contour over elements

2020-01-15 Thread David Eaton 
I understand the C0 element is piecewise linear across elements. However, I did 
not experience the same issue in my own C++ code while I use C0 element with 
the Petrov Galerkin stabilization terms. Actually, I am very confused at this 
point. How could I get rid of it while using C0 element?

Thanks
D.


From: dealii@googlegroups.com  on behalf of Wolfgang 
Bangerth 
Sent: Thursday, January 16, 2020, 12:05 AM
To: dealii@googlegroups.com
Subject: Re: [deal.II] discontinous contour over elements

On 1/14/20 10:04 PM, David Eaton wrote:
>
> Thank you for your suggestions. I am going to take a look at Lethe and compare
> with my implementation. In stabilized formulation, I used quadrilateral
> element, instead of P2 P1 Taylor-Hood element. The used element is only C0
> element. I also did not expect such a discontinuity between elements. Although
> I use P2 P1 Taylor-Hood element without stabilization terms, the discontinuity
> is still there. Probably I made mistakes somewhere  in setup. I also suspect
> that my solution is not converged. After taking a small relative tolerance
> 1e-8, the discontinuity still appears.

David, you did not understand what we were saying: If you use C0 elements
(think, piecewise linear) and you take derivatives to compute the vorticity,
then you automatically get a discontinuous function. That has nothing to do
with stabilization, solver tolerances, etc. It's just a consequence of the
fact that C0 elements and their shape functions have kinks and consequently
their derivatives are discontinuous.

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/f03c64fc-3736-340a-1cdc-f24749fb413f%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/SL2PR01MB2553E34853B991625FBB2DADC3370%40SL2PR01MB2553.apcprd01.prod.exchangelabs.com.

Re: [deal.II] discontinous contour over elements

2020-01-15 Thread Wolfgang Bangerth 
On 1/14/20 10:04 PM, David Eaton wrote:
> 
> Thank you for your suggestions. I am going to take a look at Lethe and 
> compare 
> with my implementation. In stabilized formulation, I used quadrilateral 
> element, instead of P2 P1 Taylor-Hood element. The used element is only C0 
> element. I also did not expect such a discontinuity between elements. 
> Although 
> I use P2 P1 Taylor-Hood element without stabilization terms, the 
> discontinuity 
> is still there. Probably I made mistakes somewhere  in setup. I also suspect 
> that my solution is not converged. After taking a small relative tolerance 
> 1e-8, the discontinuity still appears.

David, you did not understand what we were saying: If you use C0 elements 
(think, piecewise linear) and you take derivatives to compute the vorticity, 
then you automatically get a discontinuous function. That has nothing to do 
with stabilization, solver tolerances, etc. It's just a consequence of the 
fact that C0 elements and their shape functions have kinks and consequently 
their derivatives are discontinuous.

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/f03c64fc-3736-340a-1cdc-f24749fb413f%40colostate.edu.

[deal.II] Different shape representations with manifolds on the same triangulation

2020-01-15 Thread Juan Carlos Araujo Cabarcas 
Dear all, 

I would like your guidance on how to perform the assembly of different
shape representations on the same triangulation and on the same loop 
through cells.

Let me try to explain a bit more.

I have designed a grid containing two concentric squares (or circles).
Additionally, I have a parametric representation T(alpha) that maps points 
(x,y) on the edges 
of the inner square (or circle) to points (x',y') on the edges of a desired 
shape 
contained inside the outer square (or circle). This is 
(x',y')=T(alpha;x,y), where 
the representation T depends on a parameter alpha.

Now assume we have N different shapes that are labelled with the index m.
We then have N parameters alpha_m, with m=1,2,...N.
Then in the loop through the triangulation cells in the assembly process, I 
would like
to be able for each cell to loop through different alpha_m, in order to 
generate the 
local FE matrices and load vectors corresponding to each of the 
modified triangulations corresponding to each of the N shapes.


My guess is that somehow I should modify the following line:

const FEValues _v = hp_fe_v.get_present_fe_values ();


and pass instead a FEValues that has been generated with the desired 
alpha_m.

Below I provide a more complete sketch on how my code looks like.

I am quite hesitant on where to start, and I would appreciate your insights 
on 
how to achieve having an assemly routine where I can pass different alpha_m 
per each 
cell in the loop through cells.
I am grateful for any advise or hint on how to achieve this.

Thanks in advance,
  Juan Carlos Araújo, PhD



The way I work with one shape is the following:

// Deaclare my environment
PolarShapeManifold manifold;
PolarManifold polar;
TransfiniteInterpolationManifold inner_manifold;

Triangulation triangulation;

hp::DoFHandlerdof_handler;
hp::FECollectionfe_collection;
hp::MappingCollection mapping_collection;

PETScWrappers::SparseMatrixsystem_matrix;
PETScWrappers::MPI::Vector solution;

// Constructor
manifold ( alpha), dof_handler (triangulation),...

// setup_system
typename hp::DoFHandler::active_cell_iterator
cell = dof_handler.begin_active(),
endc = dof_handler.end();
for (; cell!=endc; ++cell) { 
cell->set_active_fe_index ( ... );
}
dof_handler.distribute_dofs (fe_collection);
...

// create_coarse_grid 
concentric_disks_inner_shape ( triangulation );
tria.set_manifold ( manifold_label, manifold ); 
tria.set_manifold ( layer_label, polar );

unsigned int not_curved = 1;
inner_manifold.initialize(triangulation);
tria.set_manifold (not_curved, inner_manifold);
...

// Assembly

hp::FEValues hp_fe_v ( mapping_collection,

fe_collection,
quadrature_collection,
update_values|  update_gradients |
update_quadrature_points  | 
 update_JxW_values);

FullMatrix   cell_system;
Vectorcell_rhs;

std::vector local_dof_indices;

typename hp::DoFHandler::active_cell_iterator
cell = dof_handler.begin_active(),
endc = dof_handler.end();
for (; cell!=endc; ++cell) {
const unsigned int   dofs_per_cell = cell->get_fe().dofs_per_cell;

cell_system.reinit (dofs_per_cell, dofs_per_cell);
cell_rhs.reinit (dofs_per_cell);
...

hp_fe_v.reinit (cell);
const FEValues _v = hp_fe_v.get_present_fe_values ();

for (unsigned int q_point=0; q_point x = fe_v.quadrature_point (q_point);
  
  for (unsigned int i=0; ihttp://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/92fc053d-ca4e-4830-bca0-2f211a9ba384%40googlegroups.com.