On 16.06.2011 16:07, Matthew Knepley wrote: > On Thu, Jun 16, 2011 at 1:43 PM, Alexander Grayver > <agrayver at gfz-potsdam.de <mailto:agrayver at gfz-potsdam.de>> wrote: > > >> 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? :) > > > Do you ever actually use A? If not, why not just build the transformed > operator directly?
Two reasons: 1. It's wrong from the architectural point view. At the point where I compute A I prefer know nothing about transformation, just working with original physical parameters is much more natural. 2. This original A might be used in the future. I really suggested that would be easier to apply this transformation to the whole matrix. If it's impossible than I have to find another way. Just want to be sure of that. > > Matt > > 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. > > > > > -- > What most experimenters take for granted before they begin their > experiments is infinitely more interesting than any results to which > their experiments lead. > -- Norbert Wiener -------------- next part -------------- An HTML attachment was scrubbed... URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20110616/68a8bba8/attachment.htm>
