[petsc-users] Getting access to matrix rows without and setting values simultaneously

[Alexander Grayver]

 >> 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>