I would like to add a cross-gradient term in the joint inversion and need to calculate the spatial difference of the model, but I am not so familiar with the operators. Specifically, for a 2D model, the tensor mesh contains cellGradx and cellGrady to calculate the cell centered gradient, but it takes us to cell faces. What if I need to compute the spatial difference? For example, for 1D vector x=[0,1,2,3,4], the spatial difference could be D.dot(x)=[1,1,1,1] if the dimension of D is 4x5. But I think typically, we will define D as 5x5 so we will still get a 5x1 vector. I didn’t find a differential operator that can achieve this. Could you please help to figure it out? Thanks.
This does indeed means that when taking the gradient in discretize, the values go from the cell-centers to the faces of the mesh. discretize has many operator that you can see by taping “Tab”, like mesh.[Tab]. Most will tell you in the name what they take (cell-centers: CC), and where they go (Faces or Nodes).
Also, a research group in Houston with Professor Jiajia Sun is currently working on bringing an implementation of the cross-gradient to SimPEG (see Pull Request 943). It might be worth for you to contact them directly in person (maybe through the SimPEG Slack group).
Related more specifically to the desired size of your operator, I am not sure of what is the issue you encounter. If the gradients of the two models have the same size, then it should not cause an issue to compute the cross-gradient, even if the gradient have a different size than the original model.
If this is a necessity, you could try to first take the gradient from CC to F, and then average it from F to CC? (note that this is not a recommendation from my side to use such operator, simply an example on how to build it): mesh.aveF2CC * (mesh.cellGrad * model)
Thank you for your reply. The issue is not about the size of the two models. They are the same. It’s because I need to compute the gradient along different dimensions. if nCx is not equal to nCy, the dimension of mesh.cellGradx * model will not equal to that of mesh.cellGrady * model because the number of faces along x and y directions is different. Actually, I tried to use mesh.aveFx2CC * (mesh.cellGradx * model) to compute the spatial gradient along x direction, but when I go to the detail of the “equivalent” differential operator Dx=mesh.aveFx2CC * mesh.cellGradx, it is different from what a general differential matrix (with non-zero adjacent -1 and 1 in each row) looks like. As marked with the red lines in the following picture, there are zeros between the non-zero elements. It seems that Dr. Sun’s group adopted a similar way to compute Dx or Dy, but I am not sure what they finally got. I think it looks a little bit unreasonable to use a differential matrix like this. So, overall this is the way on SimPEG about how to start from CC and get back to CC, right? Even though there is some difference from the general form of the differential matrix.
Ok, I understand now what your issue is with the differences of sizes between the gradients in the X and Y direction, thanks for the clarification.
As for if Dx is a reasonable operator, I can indeed see issues arising as its null space includes some periodic signals. For example, if you have a periodic signal like m = [-1, 1,-1, 1, [...] ], then Dx * m = [0, 0, 0 [..]]. So you get a 0 gradient while it is clearly not true. How this aspect is solved in practice, I simply do not know, it is not something I have personally investigated. There should be something in the literature I hope…
Thank you so much for your help. Will go through the literature you mentioned and get back to you if I find the answer in case you need it in the future
Thank you for your information. I “extracted” the cell center to cell center differential operator (differential operator in general) from mesh.cellGrad by deleting some rows for the 2D case. For the 3D case, it’s probably not that straightforward and what you shared I think can help a lot. Thanks again for sharing.
