[petsc-users] Question on matrix assembly

[Sajid Ali]

@Matthew Knepley <knepley at gmail.com> : Thanks for the explanation on
preallocation.

>However, why not have a flag so that on the first pass you do not compute
entries, just the indices?

The matrix computes the projection of an image onto a detector so
generating this involves computing all possible ray-rectangle intersections
and computing the values only differs from computing the indices by a call
to calculate intersection lengths. The process to set up the geometry and
check for intersections is the same to generate indices and values.

So, in this case the tradeoff would be to either compute everything twice
and save on storage cost or compute everything once and use more memory
(essentially compute the matrix rows on each rank, preallocate and then set
the matrix values).

@Stefano Zampini <stefano.zampini at gmail.com> : Yes, I only need MatMult and
MatMultTranspose in the TAO objective/gradient evaluation but in the
current state it's cheaper to use a matrix instead of computing the
intersections for each objective/gradient evaluation. About ~70% of the
application time is spent in MatMult and MatMultTranspose so we're hoping
that this would benefit from running on GPU's.

Thanks for the pointer to MatShell, implementing a matrix free method is
something we might pursue in the future.

-- 
Sajid Ali | PhD Candidate
Applied Physics
Northwestern University
s-sajid-ali.github.io
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20200814/f1a0022f/attachment.html>