Rotated Gradients #1167
Rotated Gradients #1167
Conversation
Also adds more documentation with examples
replace with n_cells
|
Ah, forgot there were changes to the DC simulation on here as well, I'm going to parse those out and put them on a separate PR |
|
|
|
A comparison of an inversion run with the current smoothness operator on a tree mesh, and this smoothness operator on a triangular mesh. |
Codecov ReportAttention: Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## main #1167 +/- ##
==========================================
+ Coverage 82.32% 82.36% +0.04%
==========================================
Files 165 166 +1
Lines 25263 25379 +116
==========================================
+ Hits 20798 20904 +106
- Misses 4465 4475 +10 ☔ View full report in Codecov by Sentry. |
|
Just having a thought here... Is it actually necessary to provide orthogonal vectors for the regularization directions, If someone wanted to regularize in two directions that were not orthogonal, I think this formulation still holds, because of the way the anisotropic weighting parameter is constructed. |
…into anisotropy_inversion_support
SimPEG/regularization/rotated.py
Outdated
| G = self.cell_gradient | ||
| M_f = self.W | ||
| r = G @ (self.mapping * (self._delta_m(m))) | ||
| return 0.5 * r @ M_f @ r |
There was a problem hiding this comment.
Here we only have one M_f with M_f = self.W, don't we miss another to have W^TW ?
There was a problem hiding this comment.
(same comment for deriv and deriv2)
There was a problem hiding this comment.
Actually no, the way it works here is that M_f=W.T@W in the standard way of describing the regularization. Thus since I use discretize to construct the anisotropic inner product matrix. So I just set W=M_f since there’s a lot of machinery already hooked into updated the weighting matrix on things like weighting function changes.
There was a problem hiding this comment.
Thanks. I do admit I find this extremely confusing and important to document here.
There was a problem hiding this comment.
M_f = self.W; then should be M_f = self.WTW for clarity
SimPEG/regularization/rotated.py
Outdated
| """ | ||
| if getattr(self, "_W", None) is None: | ||
| mesh = self.regularization_mesh.mesh | ||
| cell_weights = np.ones(len(mesh)) |
There was a problem hiding this comment.
@jcapriot : I don't think this works with active_cells because len(mesh) here will be the size of the total mesh (active + inactive), but later, you multiply those cell_weights by values (such as sensitivity weights) which will be the shape of the number of inactive cells.
See the error I am getting:
There was a problem hiding this comment.
PS: the issue only appears if you are setting weights AND have active_cells.
|
Hey @thibaut-kobold @JKutt @johnweis0480 @lheagy Has promised, here's a VERY old notebook that I had used for plots in my thesis for sparse rotated gradients. I tried updating it, but turned out to be a lot of work.
It relies on a simpeg/SimPEG/Regularization.py Line 2100 in ac1a092 Technically you could just copy that part only and overload the gradient operators on the reg classes to try it. This stuff was never published, but someone should do something about it. Enjoy. |
SimPEG/regularization/rotated.py
Outdated
| n_faces = mesh.n_faces | ||
| n_cells = self.regularization_mesh.n_cells | ||
| cell_weights = np.ones(n_cells) | ||
| face_weights = np.ones(n_faces) |
There was a problem hiding this comment.
@jcapriot : curious if those face weights will allow for sparsity? For what I gather, it won't, correct?
There was a problem hiding this comment.
They'd likely allow for it to be sparse in x-y-z orientations, but I don't think they would work properly for making it sparse in a rotated direction. It might be approximately that on minimization but I'm not sure on that front either because I haven't tried it.
Constructs a (optionally) rotated gradient operator on general meshes.
This makes use of discretize's
get_face_inner_productoperator to construct the regularization as an inner product of cell gradients. Essentially the weighting parameters are an anisotropic model used to combine the gradients.This operator also works on general mesh types to construct the smoothness regularization operators on unstructured meshes as well.
Things to check (in addition to the implementation).
sqrt-ed and applied on both sides of the current_Wmatrix.)totaloption is handled, but it will require some experimentation about how this applies to rotated systems.