[petsc-users] Getting access to matrix rows without and setting values simultaneously
Alexander Grayver
agrayver at gfz-potsdam.de
Thu Jun 16 08:43:46 CDT 2011
>> Could you explain more about the task you're trying to do.
Well, I can try. That is pretty specific task, I wouldn't like to go
deep into details.
Let's say we have a discretized model of some physical parameter m (say
acoustic velocity). Number of model parameters is N. We need to take M
measurements d within model (say time traveling of acoustic wave) based
on M different sets of receiver/source positions.
In our case N >> M (e.g. 10^7 >> 10^3).
We have operator F (which might linear or not) that defines relationship
between m and d:
d=F(m)
The operator F is just a system of equations actually. I have no problem
now to solve it for any m.
What I need now is to compute the variation of this operator:
A = \frac{\partial d_i}{\partial m_j}, for all i=1..M, j=1..N
After some maths this calculation could be reformulated as a triple product:
A_i = C*F(m)^-1*v
Where C some sparse matrix.
Once the A is computed I want to solve another problem:
(A'*A)b=A'*r
Which is a Gauss-Newton system now,
A'*A is truncated Hessian,
A'*r is gradient,
r = d - d_obs, where d_obs is the real observed data.
b -- model change which have to applied to original model m in order to
explain your observed data better
The latter problem I want to solve using petsc matrix-free formulation
and some iterative solver (haven't decided yet which, could you advice
one?).
But the POINT here is that the latter problem must be solved not in the
original m-space, but in transformed x-space. For that we need another A:
A_tr = \frac{\partial d_i}{\partial x_j}
However, using chain rule you can represent A_tr in terms of A:
A = \frac{\partial d_i}{\partial m_j} * \frac{\partial m_j}{\partial x_j}
\frac{\partial m_j}{\partial x_j} - is the exactly transformation I want
to apply to matrix A. It's a simple scalar expression, but has to be
applied to each element of A.
Did it help? :)
Regards,
Alexander
On 16.06.2011 14:55, Jed Brown wrote:
> On Thu, Jun 16, 2011 at 14:48, Alexander Grayver
> <agrayver at gfz-potsdam.de <mailto:agrayver at gfz-potsdam.de>> wrote:
>
> by the grid I meant that this matrix is 2d array, not a real grid
> of some physical parameters or whatever and it's also not a linear
> operator itself.
>
>
> Could you explain more about the task you're trying to do. Pseudocode
> or Matlab for the whole process would be useful. There might be a
> better way to do this without using Mat to store something that is not
> a linear operator.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20110616/478fc7c3/attachment.htm>
More information about the petsc-users
mailing list