[aspect-devel] Advice for optimizing particle advection for 2D Spherical Shell Models

[Rene Gassmoeller]

Hi Harsha,

I would dispute that 'scaling' is the right word for this problem. Since
your solution time scales linearly with the number of particles this is
an optimal scaling (since each particle produces the same effort). But
you are right, computing the velocity of each particle (and in
particular computing the location in the reference cell for higher order
mappings / curved geometries) is a major part of the whole particle cost.

Unfortunately, I do not think there is an easy way around that, you
either need to compute the particle's location in the reference cell, or
need to project the real velocity into the coordinate system of the
reference cell (which is at least as expensive, and introduces another
source of error). You are right that transforming the velocity at the
reference points might save you some time, but I am not yet convinced
that this would also save time for higher order integration schemes like
RK2 where you would need to project the old and the new velocities.
Additionally, you would still need to recompute the reference location
when the particle moves to another cell. The other options to speed up
this process are equally unattractive: You could use a first order
mapping for the cells to make the inversion faster (and vastly decrease
the accuracy of the Stokes solution near the boundaries), you could use
Q1P0 elements that at least make the evaluation of the velocity faster
once you have the position in the reference cell (thereby picking up a
lot of velocity divergence at the particle positions), or you could
rewrite the organization and whole particle algorithm in a
mesh-independent way, which is a major effort. During writing that
algorithm I spend quite some time thinking about the same problems, but
did not come up with a better solution so far. If you find one, feel
free to let me know, I would appreciate a speed-up in that place a lot ;-)

I annotated some of your questions below to the best of my knowledge.

Best,

Rene


On 09/06/2017 01:54 AM, Harsha Lokavarapu wrote:
>  
> Anyways, as a starting point, we decided to pick the simplest case,
> namely the Euler integration scheme to test our idea. Currently the
> physical location of the particle is updated using Euler. Similarly,
> we would like to introduce a second Euler step to update the reference
> location of the particle. However, this requires the velocity solution
> at the support points in the reference cell. If we are using a Q2
> velocity element, we can get away with 9 transformations.

True, but what about higher order integration schemes?

>
> After, burrowing through the Deal.II documentation, I have found
> implemented functionality for transforming vectors from unit to real
> but not the other way around. Is there a black box solution that
> already exists within Deal.II which converts vectors from real to unit? 
>
Sorry, I do not know one, but there still might be one.

> We also looked through the implementation of
> transform_from_real_to_unit_cell() which uses a newton iteration to
> transform a point in physical space to a point in the unit space. Will
> a similar technique be necessary to transform a vector from real to unit? 
>
I would guess so.

> Also, are there any locations within Aspect that involves working with
> the velocity in the reference cell?
>
No, so far there are no locations I am aware of.

> For a box domain, the transform_from_real_to_unit is orders faster
> than in the spherical shell domain. Is this because we have a better
> initial guess to the newton iteration?
>
Probably, but also because the nonlinearity is much less strong. As far
as I remember for an undeformed box we do not even need a Newton
iteration, because an analytical expression is available for
MappingEulerian (after all it is just a scaling and movement of the
reference cell, no deformation).

> We (I) have been stumped on this for close to a week now. Any leads
> are highly appreciated. 
>

I have spent a significant amount of time on the same problems two years
ago, I also appreciate any leads :-)

> Thank You,
> Harsha
>
>
> _______________________________________________
> Aspect-devel mailing list
> Aspect-devel at geodynamics.org
> http://lists.geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel

-- 
Rene Gassmoeller
http://www.math.colostate.edu/~gassmoel/

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.geodynamics.org/pipermail/aspect-devel/attachments/20170906/d7140a97/attachment-0001.html>